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Abstract 


In  this  paper  we  consider  the  approximate  solution  of  a  cl2Lss  of  second 
order  elliptic  equations  with  rough  coefficients.  Problems  of  the  type 
considered  arise  in  the  analysis  of  unidirectional  composites,  where  the 
coefficients  represent  the  properties  of  the  material.  We  present  several 
methods  for  this  class  of  problems,  and  show  that  they  have  the  same  accuracy 
ais  usual  methods  have  for  problems  with  smooth  coefficients.  We  refer  to  the 
methods  as  special  finite  elements  methods  because  they  are  of  finite  element 
type  but  employ  special  shape  functions,  chosen  to  effectively  model  the 
unknown  solution. 
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Special  Finite  Element  Methods  for  a  Class  of 


Second  Order  Elliptic  Problems  with  Rough  Coefficients 

1.  Introduction 

In  this  paper  we  consider  the  approximate  solution  of  a  class  of  second 
order,  two  dimensional  elliptic  boundary  value  problems  with  rough  or  highly 
oscillating  coefficients.  We  apply  an  approach  proposed  by  BabuSka  and  Osborn 
[5]  for  the  approximate  solution  of  problems  with  rough  input  data.  This 
approach  was  applied  to  one  dimensional  problems  in  BabuSka  and  Osborn  [4]. 
Specifically,  we  consider  boundary  value  problems  of  the  form 

(1.  1) 

■Lu(x,y)  3  -|^(a(x,y)|^(x,y)) -|y(a(x,y)|yu(x,y))  =  f(x,y)  V  (x,y)  e  n 
,u(x,y)  =0  V  (x,y)  €  30, 

2  2 

where  0  is  a  bounded  domain  in  IR  ,  f  is  a  function  in  L  (£1),  auid  where 
the  function  a  €  L”(n)  satisfies 

(1.2)  0  <  a  ^  a(x,y)  S  p  <  a>  V  (x,y)  e  0, 

where  a  and  P  are  constants.  Throughout  most  of  the  paper  we  will  also 
assume  that  a(x,y)  locally  varies  sharply  in  at  most  one  direction,  a 
requirement  on  the  coefficient  a  that  will  be  made  precise  later  (see 
Remarks  2.1  and  4.1);  such  coefficients  will  also  be  called  (curvilineair  or 
straight  line)  unidirectional.  If  the  coefficient  a  is  rough,  then  the 
solution  u  to  (1.1)  will  also  be  rough;  to  be  specific,  u  will  not  in 
general  be  in  H^(n)  and  may  not  be  in  H^**^(n)  for  euiy  c  >  0. 

Problems  of  this  type  arise  in  many  applications;  we  will  be  especially 
concerned  with  applications  to  unidirectional  composite  materials  (briefly, 
composites).  In  these  applications  the  coefficient  a(x,y)  represents  the 


prop)erties  of  the  material,  and  chainges  abruptly.  We  will  also  be  interested 
in  problems  in  which  a(x,y)  changes  smoothly  but  rapidly.  We  will  take  the 
liberty  of  referring  to  both  types  of  problems  as  composites.  In  Figs. 

1. 1-1.4  we  show  some  typical  configurations  for  unidirectional  composites. 

In  these  figures  the  coefficient  is  constant  or  is  cheinging  slowly  along  the 
lines  or  curves  and  is  changing  sharply  in  the  transverse  direction;  the 


absence  of  lines  in  a  portion  of  the  material  indicates  a  constsuit  or  a  slowly 
varying  coefficient.  We  can  interpret  the  lines  2is  fibers  in  the  composite. 
This  interpretation  is,  of  course,  symbolic  for  problems  in  which  a(x,y) 
changes  smoothly  but  rapidly. 

Fig.  1.1  shows  a  straight  line  unidirectional  composite  or  coefficient 
and  Fig.  1.2  shows  the  cross-section  of  a  tubular  composite.  Figs.  1.3a, b 
show  reinforced  panels.  The  area  A  in  Fig.  1.3b  indicates  a  region  in  which 


Are*  A 


Fig.  1.3a.  Reinforced  Panel  Fig.  1.3b.  Reinforced  Painel 

a(x,y)  is  smoothly  becoming  a  constant  or  a  smooth  function.  We  refer  to  the 
materials  or  the  coefficients  in  Figs.  1.2-1. 3  ais  curvilinear  unidirectional. 
Still  more  complicated  problems  can  be  considered.  An  example  is  shown  in 
Fig.  1.4.  While  Figs.  1. 1-1.3  depict  problems  in  which  a(x,y)  is  unidirec¬ 
tional  everywhere  in  the  domain,  the  problem  depicted  in  Fig.  1.4  corresponds 
to  a  coefficient  that  is  unidirectional  except  on  certain  lines.  Hence  we 


Fig.  1.4  Irregular  Unidirectional 
Composite 

call  it  an  irregular  unidirectional  composite.  This  type  of  problem  will  be 
addressed  in  detail  in  a  forthcoming  paper.  We  note  that  certain  interface 
problems  car  be  naturally  treated  as  problems  of  composites.  With  this 
approach  it  is  not  necessary  to  fit  the  interface  with  the  finite  elements,  as 
is  done  with  the  standard  approach. 

A  finite  element  method  is  obtained  by  restricting  the  weak  formulation 
of  problem  (1.1), 

ueHgCn) 

(1.3)  •  , 

B(u,v)  3  a  grad  u  •  grad  v  dxdy  =  fv  dxdy  V  v  e  H  (£1), 

I-  '’n  "  ''n 

to  finite  dimensional  trial  and  test  spaces.  The  outline  of  the  approach 
given  in  [5]  is  ais  follows: 

1)  Characterize  the  space  of  solutions  corresponding  to  the  space  of 

2 

right-hand  sides  (in  our  case  we  suppose  f  e  L  (£1)).  This  will  involve  a 
regulsu'ity  result.  Although  regularity  results  are  well-known  for  elliptic 
problems  with  smooth  coefficients,  they  are  not  available  in  a  direct  form  for 
our  problem.  Such  results  will  be  discussed  in  Section  2. 
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2)  Select  trial  spaces  that  have  good  approximation  properties.  The 
approximation  properties  of  the  trial  functions  or  shape  functions  are 
directly  tied  to  the  regularity  of  the  solution.  For  example,  if  the  solution 
u  of  (1.1)  is  not  in  H^(n),  then  it  is  well-knovm  that  the  usual  finite 
element  method  based  on  piecewise  llneeu'  approximating  functions  produces 
inaccurate  results.  The  problem  of  selecting  optimal  trial  functions  is  not 
simple;  in  practice,  one  would  like  to  find  a  trial  space  that  performs  almost 
as  well  as  the  optimal  one  but  that  can  be  reasonably  implemented.  We  use 
the  phrase  special  finite  element  methods  to  denote  methods  with  this  type  of 
special  shape  functions. 

3)  Select  a  test  space  so  as  to  ensure  the  inf -sup  (or  stability)  con¬ 
dition  is  satisfied  and  so  that  the  resulting  finite  element  method  can  be 
reasonably  implemented. 

We  will  use  this  approach  to  design  methods  of  finite  element  type  which 
will  yield,  roughly  speaking,  the  same  accuracy  as  the  usual  finite  element 
method  when  a  is  smooth,  but  strikingly  improved  accuracy  when  a  is  rough. 

The  orgzuiization  of  the  paper  is  ais  follows.  In  Section  2  we  present  the 
regularity  results  needed  for  the  problems  we  eire  dealing  with.  (Although 
Theorem  2.4  presents  a  regularity  result  for  problems  of  the  type  depicted  in 
Fig.  1.4,  the  complete  treatment  of  which  will,  as  mentioned  above,  be  the 
subject  of  a  forthcoming  paper,  we  have  included  it  for  completeness. )  Then 
we  will  propose  and  analyze  several  methods  to  solve  problem  (1.1)  in  the 
special  case  in  which  0  =  11^=  (0,1)  x  (0,1)  and  a(x,y)  =  a(x)  is  a 
function  of  x  only.  This  study  is  carried  out  in  Section  3,  where  we 
propose  three  distinct  approximation  methods.  A  function  a(x,y)  =  a(x)  of 
X  only  is  an  example  of  a  function  that  locally  varies  sharply  in  at  most  one 
direction;  in  fact,  such  an  a(x,y)  globally  vairies  sharply  in  at  most  one 
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direction.  a(x,y)  cam  also  be  referred  to  as  straight  line  unidirectional 

(see  Fig.  1.1).  In  Section  4  we  present  a  further  development  of  two  of  the 

methods  from  Section  3  in  order  to  treat  problems  of  the  type  depicted  in 

Figs.  1.2  and  1.3  with  curvilinear  unidirectional  coefficients. 

As  noted  above  the  approach  presented  in  this  paper  is  thoroughly  studied 

in  the  one-dimensional  case  in  [4].  Techniques  similar  to  special  elements 

were  used  in  Ciarlet,  Natterer,  and  Varga  [8]  amd  in  Crouzeix  and  Thomas  [9] 

to  handle  degenerate  one-dimensional  elliptic  problems.  We  also  mention  the 

recent  work  of  Moussaoui  amd  Ziami  [16],  which  deals  with  the  same  kind  of 

problems  with  a  method  similar  to  our  Method  I,  presented  in  Subsection  3.1. 

Finally  we  mention  the  papers  [3,  17,  18],  which  are  related  to  our  approach. 

Throughout  the  paper,  we  will  use  the  L^(n)-bamed  Sobolev  spaces  H^(n), 

consisting  of  functions  with  partial  derivatives  of  order  less  than  or  equal 
2 

to  k  in  L  (n).  These  spaces  are  equipped  with  the  norms  and  semi-norms 


We  will  also  use  the  spaces  H^(n)  for  fractional  k.  consists  of 

those  functions  in  H^(n)  which  vanish  on  30.  We  will  also  use  the  space 
H  ^(£1)  =  [Hq(£1)]'.  Throughout,  C  will  denote  a  generic  constant.  When  we 
say  "there  exists  C  =  C(a,8),"  we  mean  that  C  depends  on  the  coefficient 
a(x,y)  only  through  its  upper  amd  lower  bounds  a  and  8  (cf.  (1.2)). 

The  authors  would  like  to  thank  Professor  L.  C.  Evans  for  calling  the 
theorem  of  S.  N.  Berstein  used  in  Section  2  to  their  attention. 
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2.  Regularity  Results 


It  is  clear  that  problem  (1.1)  has  a  unique  (weak)  solution  in  Hq(£2); 
cf.  (1.3).  This  is  an  inunediate  consequence  of  the  Laoc-Milgram  theorem. 
Furthermore, 


|u||i^^.C(a)||fllo^^. 


1 

But  if  a(x,y)  is  rough,  then  u  may  not  be  in  H  (fl)  for  auiy  e  >  0,  and 
therefore  we  cannot  expect  any  reasonable  rate  of  convergence  for  the  usual 
finite  element  method.  Nevertheless,  as  a  consequence  of  the  assumption  that 
a(x,y)  is  unidirectional,  the  solution  u  does  satisfy  a  regularity  property 
that  can  be  employed  in  the  derivation  of  sui  approximation  method  for  (1.1) 
with  a  good  rate  of  convergence,  even  though  a(x,y)  is  rough. 

It  is  the  purpose  of  this  section  to  prove  such  regulau~ity  results,  first 
for  the  model  problem  consisting  of  (1.1)  with  £1  =  £1^  =  (0,1)  x  (0,1)  and 
the  coefficient  a(x,y)  satisfying  a(x,y)  =  a(x),  i.e.,  with  a  straight  line 
unidirectional  coefficient  (cf.  Fig.  1.1),  and  then  for  the  more  general 
problem  with  a  curvilinear  unidirectional  coefficient  (cf.  Figs  1.2-1. 3). 
Finally  we  prove  a  regularity  result  that  will  be  applied  in  a  forthcoming 
paper  to  a  problem  with  an  irregular  unidirectional  coefficient  (cf.  Fig. 

1.4).  Our  main  tool  is  a  theorem  of  Bernstein  [6],  [13,  Section  3.17]  for 
elliptic  equations  in  non-divergence  form,  which  we  now  state. 

Consider  the  problem 


(2.  1) 


2  2  2 
a  u  a  u  a  u 

"  ~  12axay  "  ^22gy2 

u  =  0  on  a£l, 

V.  ’ 


f  in  £1 


2 

where  £1  is  a  bounded  convex  domain  in  IR 
2 

C  boundary  a£l  and  where  the  functions  a 


with  a  Lipschitz  and  piecewise 


ij 


€  L  (£1)  satisfy 
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(2.2)  y 


2  2  2 

i=l  i.J=l  i=l 

with  agj  =  aj2.  where  v  and  (i  are  positive  constamts.  Note  that  the 
equation  in  (2.1)  is  in  non-divergence  form. 

2 

Theorem  2.  1  (Bernstein).  For  each  f  6  L  (£1),  problem  (2.1)  has  a  unique 
solution  u  €  H^(£l)  A  Hq(£1).  Furthermore,  there  is  a  constant  C(y,p), 
depending  on  v  and  p  but  independent  of  f,  such  that 

(2.3)  llullj  „  £  C(,..M)l|f|lo  s,. 

Our  hypothesis  on  £1  is  not  identical  to  the  one  in  [13].  To  prove 
that  (2.3)  is  still  valid  for  such  a  domain,  one  can  use  the  a  priori 
estimates  given  in  [11,  Section  3.1]. 

The  first  application  of  Bernstein’s  Theorem  will  give  a  regularity 
result  for  problem  (1.1)  when  £1  =  £1^  and  a(x,y)  =  a(x).  Corresponding  to 
problem  (1.1),  with  this  assumption,  we  define  the  space 


2 

Theorem  2.2.  Suppose  £1  =  £1^  and  a(x,y)  =  a(x).  Then  for  each  f  e  L  (Q) 
the  solution  u  of  (1.1)  is  in  ^^(£1)  n  H^(£l).  Furthermore,  th<  re  is  a 
constauit  C  =  C(a,p),  depending  on  a  and  p  but  independent  of  f,  such 
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that 


(2.6) 


|ullL.n  " 


Proof.  Let  u  be  the  \inique  solution  to  (1.1)  in  Hq(£1).  We  introduce  the 
chsuige  of  variables  or  mapping 


(2.7) 

and  the  notation 

(2.8) 


x(x)  = 


ds 

J.(s) 


.  y(y)  =  y 


u(x(x),y(y))  =  u(x,y),  (x,y)  e  Q. 

.1 


The  mapping  (2.7)  maps  the  domain  Q  onto 


n  =  (0,  f 

J  r 


ds 

a(s) 


) X ( 0 , 1 ) .  We  see 


that  u  €  H^(n)  if  and  only  if  u  €  H^(n),  af^,  e  H^(n),  which  is 

OX  oy 

equivalent  to  u  €  H^(n).  We  also  note  that  the  weeik  formulation  (1.3)  of 
(1.1)  is  transformed  into 


(2.9) 


u  € 


f  fdu  dv  ~  du  3v1  r  w  -  ,,lf?.x 

- +  a - dxdy  =  f  avdxdy  V  v  6  ( n ) . 

•'’n  ^x  dx  dy  dyJ 


The  system  (2.9)  is  simply  the  variational  formulation  of 


(2. 10) 


-2  d^u 
a  — ^  =  af 

ay 


in  n 


[u  =  0  on  dQ. 


Note  that  while  the  equation  In  (1.1)  is  in  divergence  form,  the  equation  in 

(2.10)  is  in  non-divergence  (as  well  as  in  divergence)  form.  As  a  consequence 
of  Theorem  2.1,  (2.10)  is  uniquely  solvable  in  Hq(S)  n  H^(Q)  guid 

(2.11)  ^  C(a,8)||iLf||Q 

Since  (2.10)  is  uniquely  solvable  in  H^lQ),  we  conclude  that  u,  as  defined 


9 


in  (2.8),  which  satisfies  (2.9),  coincides  with  the  solution  of  (2.10)  and 


hence  lies  in  HQ(n)  n  H^(n)  and  satisfies  (2.11).  Thus  u  €  Hj(n)  n  H^(J2), 
which  is  the  first  conclusion  in  the  theorem.  If  we  cheinge  variables  in  the 
estimate  (2.11)  to  return  to  the  original  variables,  we  obtain 


u 


L,n 


2 

u  dxdy  + 


f,au,2  ^  ,au,2') ,  . 

['a;'  *  'apl 

2  2 
f  ,3  ,  au, ,2  ,a  u  ,2  i.a  u,2'l.  . 


u^adxdy  + 


a 


a 


09 


1  |au,2  ~|9U|2l  ,~,~ 

-  I— I  +a|— I  dxdy 

a  Qy  J 

flv'^  J 


ax  ay 

2~ 

axay  ay 


dxdy 


<  max0.1)||G||2^g 


<  maxO,i)C^(a,p)||5fll2  ~ 


S  0^maxO,i)C^(a,p)llfl|^  _, 
a  u,  n 


which  is  (2.6). 


Theorem  2.2  was  proved  by  maJcing  a  global  cheuige  of  variables  and  then 
applying  the  Bernstein  result.  The  global  chauige  of  variables  exists  because 
a(x,y)  globally  varies  sharply  in  one  direction:  a(x,y)  =  a(x).  We  now 
prove  a  second  result  in  which  we  assume  the  existence  of  only  a  local  change 
of  variables  (cf.  Figs.  1.2-1. 3). 

Let  E  c  n  be  open  and  assume  that  we  have  a  system  of  orthogonal 
curvilinear  coordinates  defined  on  Z.  More  precisely,  regarding  Z 

and  (C.Tj)  we  assume 

(i)  the  functions  ^,t}  are  defined  on  Z  and  are  twice  continuously 
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differentiable, 

(ii)  (C,T>)  :  I— »Z'  is  one-to-one  emd  onto, 

Irlt^ "  ° 

(iv)  grad  ^  •  grad  tj  =  0  in  Z, 

12  12 

(v)  Z  is  a  rectangle  in  ^,7),  i.e.,  Z'  =  (Cj.Cj.)  x  and 

(vi)  Z  n  a  £1  =  xS,  in  which  case  all  edges  of  Z  are  called 

interior  edges 


the  union  of  one  or  more  edges  of  Z,  in  which  case 
these  edges  aire  called  boundary  edges  amd  the 
remaining  edges  are  called  interior  edges. 

The  union  of  the  interior  edges  is  denoted  by  E.  We  suppxjse  further  that 


(2. 12) 


a(x,y)  =  a'(C)  V  (x,y)  €  Z, 


where  we  use  the  notation,  for  any  function  w  defined  in  Z, 

w' (C(x,y),7)(x,y))  =  w(x,y),  (x,y)  €  Z. 

See  Figs.  2.1a,b  for  typical  configurations. 


Fig.  2.1a 


Fig.  2.1b 
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Theorem  2.3.  Let  u  be  the  solution  in  H^Cn)  of  (1.1),  where  we  assume 
2 

that  f  e  L  (n)  and  that  a(x,y)  satisfies  eissumptions  (1.2)  and  (2.12), 
where  I,  (^,t))  satisfies  conditions  (i)  -  (vi)  above.  Let  0  c  Z  be  open  and 
satisfy  OccZ  if  Zndn  =  0  and  dO  r\  d  'Z  c  d  Q  if  Z  n  d  Q  *  a,  and 
let  0'  be  the  image  of  0  under  the  mapping  (?,  i))  (cf.  Figs.  2.1a,b). 

Then  there  is  a  constant  C  =  C(a,0,^,i),d)  depending  on  a,  $,  tj,  euid  d 
but  independent  of  f,  such  that 


(2. 13) 


+  a 


.a^u' 

'acdT) 


2  1 


t-'O' 


,a^u' 

2 

dri^ 


d^dT) 


1/2 

<  C(a,^,?,^,,d)llfl 


o.n’ 


where  d 


dist(<9,  E),  if  E  *  0 
^1,  if  E  =  0 


Proof.  Clearly  u  (more  precisely,  u|_)  satisfies 


(2. 14) 


u  €  H^(Z) 


J  a  grad  u  •  grad  v  dxdy  =  J  fv  dxdy  V  v  €  Hq(Z). 


Introducing  the  change  of  variables  (?,7))  in  (2.14)  we  get 


(2. 15) 


f  ,  fau'  av'i  .  e:, 2^ an'  av',  ^  ,2'| 


a(x,  y) 
a(C.T7) 


dCdTj 


Now  we 

(2. 16) 


introduce  a  second  ch2Lnge  of  variables, 

r? 

f  dt 

^  ^  aTtT  ’  ^ 


where  is  the  ^-coordinate  of  the  left  edge  of  Z' .  We  will  use  the 

notation 
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w'dC?."*)).  T)(C,7)))  =  w'(C.7>),  (?,T))  e  2',  for  any  ftaiction  w' 


on  Z' 


Z'  =  image  of  Z'  under  the  mapping  (^,7}). 
Applying  this  chauige  of  variables  to  (2. 15)  we  get 


where 

(2.18)  a'(?.77)  =  Igrad  a^C^.T))  =  Igrad  t,)^ 

To  apply  Theorem  2.1,  we  need  to  introduce  Dirichlet  boundary  conditions. 
From  condition  (vi)  we  know  that  any  edge  of  Z  is  either  am  interior  edge  or 
a  boundary  edge.  Then,  through  the  correspondence  determined  by  the  mappings 
(?,  7))  and  (?,  7)),  we  will  refer  to  the  interior  and  boundary  edges  of  Z' 
and  Z' .  Now  let  e  C*(Z' )  with  ^'(^,7))  =  0  for  (^,7))  near  the 
interior  edges  of  Z' .  Then  for  v'  €  H^CZ'),  ^'v'  e  Hq(Z' )  and  we  can 
replace  v'  by  ^'v'  in  (2.17)  to  get 


or 

(2. 19) 

where 
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(2.20) 


a(C.Tj)  ““i 


3^  at 


^  a^  a?) 


Writing  w  =  u'^',  from  (2.19)  we  get 


(2.21) 


w  €  hJ(Z' ) 


r  r-  aw  av^^-,2~  ^  ^Idld^i  =  f  Fv'  d^dv  V  v'  e  hJ(Z'). 

L.  ^  ^  at  at  ^  3v  3n  J  L.  ° 


3ri  3ri 


w  €  Hq(Z' )  since  =  0  near  the  interior  edges  of  Z' .  Since  u'  is  in 

1  -  2  ~  -  i 

H  (Z'),  the  functions  are  C  ,  and  a'  =  a'(t).  we  see  that  F  is  in 

2  ~ 

L  (Z').  The  system  (2.21)  is  simply  the  variational  formulation  of 


(2.22) 


at*- ^  atJ  a^l^aV 


w  =  0  on  aZ'  . 


The  equation  in  (2.22)  can  be  formally  written  ais 


(2.23) 


-  S;?!|  -  =  F  . 

at  dv  at  at  arj  3v 


Denote  by  W  the  unique  solution  in  H^(  of 


(2.24) 


which  exists  by  Theorem  2.  1  and  which  satisfies 


(2.25) 


||W|l2^g,  ^  C(a.3.C.T,)  IIGIIq^^,. 


Now  w  €  Hq(Z')  solves  the  same  problem  formally.  We  want  to  show  that 
w  =  W,  and  hence  that  w  €  H^(Z')  and  satisfies  (2.25).  Writing  (2.24)  in 


divergence  form  we  obtain 
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a? 


(2.26) 


^  ^  at-*  dri  ^ 

-  -  a'2 

i'  - 

aw  . 

2  ^^2  aw 

'at^ 

at  at 

a?}  aT) 

.  .  ""iaw  . 

-,2^^2aw  ^^law 

~,2^Saw 

a|  a|  ai)  aij  a|  at  a^  arj 

Letting  U  =  W  -  w  and  using  (2.22)  and  (2.26)  we  see  that 


(2.27) 


at  >>  atJ  ai)  >-  aV 


^SaU  ~,2  ^^2dU  .  f,, 

=  - - a' - in  Z 

at  at  ai)  arj 

u  a  0  on  a  Z' . 

where  we  understand  the  equation  in  the  weaJc  sense.  It  will  suffice  to  show 
that  U  =  0. 

Let  T  :  H  ^(Z') — )Hq(Z')  be  the  solution  operator  corresponding  to  the 


problem  (2.22),  i.e.,  let  TF  =  w.  Then  from  (2.27)  we  get 

f 

(2.28)  U  a  T - - 

I  a? 


^"'au  .,2^-2 

-  a  - 


3  AU. 


at  at  a?}  ai 

Since  T  :  H  ^ (Z' )  — > HQ(f '  )  is  bounded  and  H*^(Z' )  is  compactly  contained 

in  H  ^(Z'),  we  see  that  A  :  Hq(Z' ) — » Hq(Z' )  is  compact.  Suppxjse  now  that 
U  5*  0.  Then  from  (2.28)  we  see  that  1  is  en  eigenvalue  of  A.  Hence  1  is 
3J\  eigenvalue  of  A  ;  let  V  be  an  associated  eigenfunction. 

We  can  choose  V'  e  H^(Z' )  n  H^Z' )  so  that  HV  -  V' ||  <  I1V||  s;,  ■ 

U  L  t  ^  1  t  b 

Then  (g' )  *  from  the  Fredholm  alternative  we  see  that  the 

problem 
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(2.29) 


(I  -  A)Z  =  V 


has  no  solution  in  Hq(2').  Recalling  the  definition  of  A 
equation  (2.29)  can  be  written 


(2.30) 


-  a'2!!i  =  V' 

^  a|  dri  aV 


we  see  that 


which  can  be  formally  written  as 


or 


^  ^  ^  _  3_f~,  avO  _  £,2a_f~,  ^1 

a|  a|  dri  an  at'-^atJ  dv^^dvJ 


a^z  ~,2-,  a^z  a  f-,  avO  ~,2a  f-,  avn 
af'"  af  ^  '  a^i"^  a^  i  an^^aj 


But  Theorem  2. 1  shows  that  this  equation  has  a  solution  Z  in  H^(Z')  n 
Hq(Z').  It  is  immediate  that  this  Z  solves  (2.30)  zind  hence  solves  (2.29), 
a  contradiction.  Therefore  we  conclude  that  U  =  0.  We  have  thus  shown  that 
w  =  W  and  hence  that  w  €  H^(Z')  and 

(2.25')  ^  C(a,B.?,T))  IICIlQ^g,. 


Now  the  function  can  be  chosen  as  a  cut-off  function  satisfying 

^'  =  1  on  5',  ID^'I  <  C(d' )  \  and  ID^^'l  £  C(d' )  where  C  is  some 
positive  constant,  zind 

^  dist(8',  E'),  if  E'  *  ^ 

[l,  if  ^ 


Then  from  the  expression  for 


IICII 


0,Z' 


G  in  (2.23)  auid  from  (2.20)  and  (2.18)  we  have 

s  c(a.B,€,Tj)(d')"^ 
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Thus  from  (2.25')  we  see  that 


(2.31)  "^'“2.8'  ^  "''"2. S'  ^  C(a.^.C.T,.d')  llfllg^jj- 

Now  changing  from  I.tj  back  to  C.  t)  (2.31)  we  obtain 


i  C(a.0,C.n.d')  llfllg  Q. 

which  yields  (2.13)  since  Cj(a,0,^,T))d'  <  d  i  C^(a,  rt)d' .  m 


Remark  2. 1.  Equation  (2.12),  with  (C.'O)  and  S  satisfying  the  conditions 
(i)-(vi),  is  basis  for  the  precise  meaning  of  the  phrase  "a  locally  varies 
sharply  in  at  most  one  direction,"  which  is  fully  formulated  in  Subsection  4.1 
(see  Remark  4.1). 


Remark  2.2.  If  the  mapping  function  ?  +  irj  is  analytic,  then  the  above 
analysis  is  simplified  since  in  this  case  the  functions  a^^  and  a^  in 
(2.18)  are  equal  to  1. 

Remark  2.3.  We  can  define  the  local  analogue  of  the  space  H^(n)  defined  in 
(2.5).  With  Z,  (?,T}),  and  0  as  in  Theorem  2.3, 


with  the  norm 


lf(0)  =  (u  :  u'  €  H^(O')  ,  a'  1^,  ^  €  H^(O')} 

OK  OJ} 


where 


,  ,2  f  f  ,,a  ,  ,au'  ,2  ,,a^u'  2  i  ,a^u',2l.,. 

'“'l.o'  h  'af'®  ar’'  **  'aws'  * i- ' rr' 

J-r,/  dri  ■' 


In  terms  of  the  semi-norm  0  ’  (2.13)  can  be  stated  as 

lul^^^  i  C(a,p,C,T,,d)|lfl|Q  jj. 
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Remark  2.4.  Theorems  2.2  and  2.3  cam  be  eaislly  generalized  to  cover 

coefficients  a(x,y)  of  the  form  a^Cxla^Cy)  amd  a^ (^laigCi)) ,  respectively, 

and  coefficients  a(x,y)  that  are  rough  in  x  but  smooth  in  y  and  a(?,7j) 

that  are  rough  in  ^  but  smooth  in  t),  respectively. 

The  method  used  in  the  proof  of  Theorem  2.3  gives  a  constant  C  which 
-2 

behaves  am  d  ,  where  d  =  dist(C?,S);  this  type  of  estimate  is  sufficient  to 

treat  problems  of  the  type  depicted  in  Figs.  1. 1-1.3.  To  analyze  the  type  of 

problem  depicted  in  Fig.  1.4  we  will  need  a  shamper  estimate.  Although 

the  complete  study  of  these  problems  is  reserved  for  a  forthcoming  paper,  we 

here  present  the  relevant  regularity  result;  we  show,  in  fact,  that  the 

-1  -2 

constant  C  behaves  as  d  Instead  of  d  by  using  a  more  refined 
analysis.  Although  the  idea  of  the  analysis  is  quite  general,  we  will  carry 
out  the  analysis  only  on  a  simple  model  problem  so  as  to  focus  on  the 
essential  idea. 

Let  0=0^=  (0,l)x(0,l)  be  the  square  with  the  boundary  dQ  composed 
of  the  four  straight  lines  F^,  i  =  1,...,4,  shown  in  Figure  2.2. 


Fig.  2.2 
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Let  VIS  further  define 


(0)  1 

HMfl)  =  {u  €  H  (fl)  :  u  =  0  on  T,  u  r„  u  T.} 

2  3  4 

and,  for  0  <  d  < 

=  {(x.y)  €  Q  :  d  <  y  <  1}. 

We  eissume  that  the  coefficient  a  €  L*(n)  is  a  function  of  x  only  and 
satisfies 

(2.32)  0  <  a  <  a(x)  <  ^  <  «  V  x  e  (0.1). 

Theorem  2.4.  For  f  e  L^(n)  and  a  6  L“(n)  satisfying  (2.32),  let  u  € 
(0) 

HMfl)  satisfy 

(2.33)  a(x)grad  u  grad  v  dxdy  =  fv  dxdy  V  v  e  11^(0). 


Then  there  exists  C  =  C(a,p),  depending  on  a  and  /3  but  independent  of  d, 
such  that 


(2.34) 


a 


(2.36)  w(x(x),  y(y))  =  w(x,y),  (x,y)  €  fl. 

For  any  function  w  in  Q,  w  is  defined  in  n  =  (0,A)x(0,l)  with  A  = 
A  tilde  means  that  we  are  in  the  domain  n. 
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Obviously  for  a  function  u  satisfying  our  hypotheses,  u,  eis  defined  by 

<0|  ^ 

(2.35)  and  (2.36),  is  in  HMfl)  and 

(2.37)  C^(a.p)  llullj^Q  i  ^ 


Furthermore  u  satisfies  the  differential  equation 

2-  2~ 

(2.38)  -  ^-^-a^(x)^-^  =  fa(x),  (x,y)  €  fi. 

ax  ay 

Here  zmd  in  the  remainder  of  this  section  we  will  interpret  such  equations  in 
the  weaik  sense.  Using  (2.37)  we  see  that  (2.34)  holds  provided  we  prove 


(2.39) 


d 


where  C  is  independent  of  d. 

The  idea  of  the  proof,  similar  to  that  of  the  proof  of  Theorem  2.3,  is 
to  apply  Bernstein’s  theorem  to  u  times  a  cut-off  function  %  chosen  so 
that  X  u  €  H|):(S).  But  to  get  a  sharper  estimate,  we  first  extract  from  u  a 


smooth  function  (0,  defined  in  (2.42))  having  almost  the  values  of  u  on  f 
This  ceuinot  be  done  in  a  simple  way,  since  we  only  know  u  e  H^(n). 

The  restrictions  of  u  to  the  edge  is  denoted  by 

w(x)  =  u(x,0),  X  €  T  =  (0,A). 


r 


Since  we  know  that 

and  ll»l, ^,0,0,1  ^  Cliailj  g, 

1/2  - 

where  11*11^/2  o  0  I  the  norm  in  the  space  (I)  (cf.  [14]),  we 

ceui  represent  w  in  terms  of  its  Fourier  series. 


w( 


k=l 
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and  the  coefficients 


(2.40) 


will  have  the  property 


k=l 


We  introduce  the  function 


(2.41) 


p(x,y)  =  sin  knx/A 


sh  kit(  l-y)/A 
sh  kn/A 


which  Is  the  solution  of 


k=l 

Ap  =  0  in  fl 
p  =  u  on  dCl  . 


Then  we  write  p  =  4>  *  where 


[d'h 


(2.42) 

II 

M 

k=l 

and 

CO 

(2.43) 

«/i(x.y)  =  y 

sh  krt  (l-y)/A 
sh  kn/A 


sin  kTtx/A 


sh  kn(l-y)/A 
sh  kn/A 


k=[d"  ]+l 


and  [d  is  the  integer  part  of  d 


(0) 


Our  method  of  proving  (2.39)  is  to  write  u  =  +  jCg*  where  e  H^(n) 

solves 


(2.44) 


-a2(x)?_|i  =  fa  in  n 

ax  ay 

Xj  =  0  on  fj. 


and  ^2  ^ 


solves 
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(2.45) 


=  0  in  n 

dx  dy 

=  0  on  fj, 


2uid  to  bound  Zj.Z2  H^(Qd)  -norm-  The  decomposition  u  =  + 

possible  since  u  satisfies  (2.38).  First  let  us  study  the  solution  of 

(2.44).  We  write  Xj  =  ^*2^,  where  e  HQ(n)  satisfies 


(2.46) 


_a^2i  _  -2  a^zi  _  ^  3^^z2  ^  ^ 

-2  ^  ~2  *^1  -2  “  ~2  ■ 

ax  ay  ax  ay 


Now  we  bound  the  L  -norm  of  F^.  Using  the  explicit  formula  (2.42)  for  <p, 
the  orthogonality  properties  of  the  functions  sin  kj  x,  coskjx,  auid  the 

A  A 

formula 


kit(l-y)  A  ^  2iik  1 

=  4ik  -  2  • 


we  obtain 


(2.47) 


I  ^—^11  ^  _  <  c  3?  k^ 

'  -2"o.n  -  ^  ^  ' 


r  ^is 

where  C  depends  on  ^  ~  J  aXTT  Independent  of  d.  Combining  (2.40) 


and  ( 2 . 47 )  we  have 


i|l"o.a^  .a- 


In  a  similaur  way  we  get 


(2.48) 


i*'2,a  ^  “'*""'i,a- 


We  also  see  that 


(2.49) 


"«i,a^'^'='i.a- 
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From  (2.46)  aJid  (2.48)  we  get 


(2.50) 


Using  Theorem  2.1,  we  see  from  (2.46)  that  e  H^(n)  amd  that 


(2.51) 


'^1**2. n  ^  c(a.p)  IIFillo^n- 


Using  (2.48)  -  (2.51)  we  see  that  is  in  H^(n)  and  satisfies 
(2.52) 


Let  us  now  consider  3;  ,  the  solution  of  (2.45).  Write  +  z  , 

^  cL  ^ 

where  z„  e  H^(f2)  satisfies 
6.  o 


(2. 53) 


a^z„  ^  ^,2. 

_ 2  _  ~2 _ 2  _  =  _  /  .  S"  'U! 

-~2  -~2  ^2  ■  _-2  --2  ’ 

3x  9v  Sx  3y 


■1, 


where  F2  €  H  (n).  Since  x^  does  not  satisfy  homogeneous  boundary 
conditions,  we  cajinot  apply  Theorem  2.1  immediately.  Let  x  €  c”(0,l) 
be  a  cut-off  function  satisfying 

j:(y)  =  1  for  d  <  y  <  1, 


and 


%(y)  =  0  for  0  <  y  <  ^, 


|;t^'^^y)  I  <  Cd  j,  y  e  (0,  1),  j  =  0,  1,2. 


The  function  x^  =  Z  Zo  «  H„(n,  .„)  is  the  solution  of 
0  £1  u  d/^ 


(2.54] 


^-a2(;)^=F3.F4  in 
3x  3y 


where  (recall  that  X2  ~ 
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(2.55) 


F_  a  — -a^(x) 

ax  3y 


(2.56) 


^■4  = 

dx  dy 


We  will  now  Ixjund  F_  and  F.  In  the  L  (fi.  ,^)-norni.  From  the  explicit 

o  4  d/d 

formula  (2.43)  for  ^  we  see  that 


'  S.  . > 


3  -7ikd/A 
e 


k=[d  ]+l 


Taking  into  account  (2.40)  and  tl^e  fact  that 


,  2  -knd/A  .  -  .-2  -ir/A  ^  .  r  ^  , 
ke  <Cd  e  Vk=(d  ]  +  !,..., 


we  obtain 


(2.57) 


2 

1 5^11  a  scd-'iiGii  g. 

dx^  ^’“d/2 


Analogously  we  get 


(2.58) 


i4|iio  a  ^  cH-‘iiSii  g, 

dy  ’d/2 


(2.59) 


ip'o.a  •  'p'o.fi  = 

3x  dy 


(2.60) 


d/2 


with  C  independent  of  d.  Using  (2.55),  (2. 57)-(2. 6C) ,  and  the  hypotheses 
on  X.  we  get 


(2.61) 


I'a'o.a,,,  ^"‘‘'“"i.a 

d/2 


Let  us  now  study  F^.  From  (2.56),  (2.45),  and  the  decomposition 
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=  yji  +  z^,  we  see  that  cam  be  written  ais 


2  2 

(2.62)  =  %(^+a^(x)^) -a^(x)(2%'(y)^+Z"(y)z5). 

ay  ay 


4  ~2 

ax 


Because  Zg  =  0  on  aS,  In  particular  in  f^,  the  Poincar6  inequality 


amserts  that 


(2.63) 


|22(x.y)|2d^dy  <  I-IZ2I1  jj. 


As  a  consequence  of  (2.53)  and  (2.59)  we  have 


(2.64) 


Now  from  (2.62),  (2.57),  (2.58),  (2.63),  and  (2.64)  we  obtain 


(2.65) 


d/2 


Theorem  2.1  can  now  be  applied  to  solution  of  problem  (2.54), 

2  ~ 

with  ^4  ^  ^  ^^d/2^  bounded  by  (2.61)  amd  (2.65),  amd  we  obtain 


(2.66) 


d  d/2 


Finally  from  u  =  z^  •••  Z2’  (2.52),  amd  (2.66)  we  get  (2.39),  which  implies 
(2.34).  I 
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3.  Special  Methods  for  Problems  with  Straight  Line  Unidirectional  Coefficients 


In  this  section  we  propose  and  einalyze  methods  based  on  special  elements 
to  solve  the  model  problem 


(3.1) 


'Lu(x,y)  3  -|^(a(x)|^(x.y) )  -  |^(a(x)|yu(x,y) )  =  f(x.y) 
u(x,y)  =0  V  (x.y)  e  dQ, 


V  (x.y)  €  n 


where  Q  =  £1^  =  (0,l)x(0,l),  f  6  L^(Q),  and  a  e  L"(n)  is  a  function  of  x 
only  and  satisfies  (1.2).  This  is  problem  (1.1)  with  a(x,y)  a  straight  line 
unidirectional  coefficient.  We  will  present  three  approximation  methods, 
prove  they  have  the  optimal  rate  of  convergence,  and  discuss  their  merits. 


3.1.  Approximation  Method  I. 

For  0  <  h  <  1,  let  be  a  triangulation  of  £1  by  (closed) 

curvilinear  triangles  T  of  diameter  i  h,  where  by  a  curvilinear  triangle 
T  c  £l  we  mean  the  pre-image  of  an  ordinary  triangle  T  c  £1  under  the  mapping 
(2.7).  Corresponding  to  we  have  a  triangulation  of  £2  by  usual 

triauigles.  We  assume  that  ^^h^0<h<l  satisfies  a  minimal  angle  condition, 

(3.2)  <<rVt€t^,  V0<h<l. 

2 

where  for  8uiy  bounded  set  S  c  R  , 


(3.3)  *^3  ~  diameter'  of  S 

and 


(3.4)  pg  =  diameter  of  the  largest  disk  contained  in  S. 

With  we  associate  the  space  of  approximating  (or  shape)  functions 

(3.5)  =  {v  €  L^(£2)  :  e  span{l,J^  i;Trr’ ^  ^  ^h’ 

V  is  continuous  at  the  nodes  of  , 

h 

V  =  0  at  the  boundary  nodes) . 


As  a  consequence  of  our  choice  for  the  ciirvi linear  triangles  T  we  see  that 

S,  c  H^Cn),  i.e.,  S,  is  conforming.  This  is  eeisily  seen  by  noting  that  the 
h  u  h 

functions  1,  f  , y  are  transformed  to  1.x, y  by  (2.7).  Consequently 

Jq  aUJ 

S.  3  (v  :  V  €  S,}  (v'  is  defined  in  (2.8)),  the  image  of  S,  under  the 
h  h  h 

mapping  (2.7),  is  the  usual  space  of  continuous  piecewise  linear  approximating 
functions  with  respect  to  and  is  conforming  since  is. 

Our  finite  element  approximation  u^  to  u  is  now  defined  by 


(3.6) 


% 


€  S. 


B(Uj^,v)  =  J  fv  dxdy  V  v  6  Sj^, 


where  B  is  defined  in  (1.3).  Uj^  is  Just  the  Ritz  approximation  to  u 
determined  by  the  variational  formulation  (1.3),  in  the  case  (3.1),  auid  the 
space  defined  in  (3.5).  Since  it  is  easily  seen  that  is  the  Ritz 

approximation  to  u  determined  by  the  variational  formulation  (2.9)  and  by 
the  space  Sj^,  we  could,  of  course,  carry  out  the  computation  and  the  analysis 
on  the  transformed  domain  Q.  We  shall  however  study  the  approximation  on  the 
original  domain  Q  since  this  approach  better  illuminates  the  more  general 
c8Lse  of  a  curvilinear  unidirectional  coefficient  studied  in  Section  4. 

It  is  immediate  that  B  is  a  bounded  bilinear  form  on  H^inixHQin). 
Furthermore,  the  stability  condition  (cf.  [1])  holds,  i.e.,  we  have 

Theorem  3.1.  There  exists  a  constant  5(a)  >  0,  independent  of  h,  such  that 
for  all  0  <  h  S  1, 


(3.7) 


inf 

veSh 


sup 

W€Sh 


|B(v.w)| 
=  1 


i  6(a) . 
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Proof.  Since  B(v,w)  is  symmetric  it  is  sufficient  to  prove  that  B  is 


coercive,  i.e.,  that 


|B(v.v)|  >  5(a)  ||vq  Vv  e  Sj^,  0  <  h  £  1. 


This  is  immediate. 


Approximability  here  involves  the  approximation  of  the  solution  u  by  a 

‘•£  iTCT 


lineeu'  combination  of  the  shape  functions  1, 


,  y  in  terms  of  which  S. 


is  defined.  Let  the  pxjints  ®  ®  ^  vertices  of  T  and  let 

P  ,P„,P„  be  the  vertices  of  T  (cf.  Fig.  3.1).  Since  the  functions 

1,  r  , y  are  transformed  to  l,x,y  by  (2.7),  we  see  that  the  inter- 

Jq  altJ 

polation  problem:  Given  numbers  find 


(3.8) 


w(x,y)  =  a  +  B 


dt  ^ 


satisfying  w(P^)  =  i  =  1,2,3,  is  uniquely  solvable. 


Suppose  u  €  H^(T).  Then  u  €  H^(t),  and  hence  u  hais  well  defined  point 
values  for  any  P  €  f.  Thus  u  has  well  defined  point  values  for  any  P  e  T, 
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and  we  define  the  span{ 1 


dt 

a(t) 


i^u  =  a  +  ^ 

J  r\ 


,  y}-interpolant  of  u  on  T  by 


We  derive  now  an  estimate  for  the  difference  u  -  d^u. 

Theorem  3.2.  There  is  a  constant  C  =  C(a,  j3),  depending  on  a,0  but 
independent  of  T  and  u,  such  that 


(3.9) 


A 


|u-d_u|  .  _  5  C;^  iul,  _  Vu€if(T). 
T  1 ,  T  L,  T 

where  h^.p^  are  defined  in  (3.3),  (3.4). 


Proof.  Using  the  transformation  (2.7),  we  have 


(3.10)  |u-d.j.uq  .j. 


|grad(u  -  d.j.u)  I'^dxdy  = 


i|£_  (u  -  d:j;u)  l^dxdy 
a  3x 


a|—  (u  -  d;ru)  l^dxdy 

f  dy 


<  maxO,^)  |u-d;ju|  j 


where  d;ju  is  the  span{  1, x, y>-interpolauit  of  u  in  the  triangle  f. 
Applying  the  usual  linear  interpolation  theorem  (cf.  [7,  p.  121]),  we  get 
the  bound 


(3. 11) 


.  .  4  ^ 

|u-d;ju|i,^  £  •^’2,7’ 


where  C  is  an  absolute  constant.  Inequality  (3.9)  is  a  consequence  of 
(3.10),  (3.11),  and  the  definition  of  the  semi-norm  |•|.  ,  with  the 

y  1 

constant  C(maLxO, -) )  ^'^^. 

a 
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We  define  now  the  Interpol  ant  of  u  e  H^(n)  by 


(3. 12) 


d,u(P)  =  u(P)  for  all  nodes  P  €  1?.  . 
■  h  h 


As  am  easy  corollary  of  Theorem  3.2,  we  can  state  our  approxlmabllity  result. 


Theorem  3.3.  There  is  a  constant  C  =  C(a,  p,<r),  depending  on  a,  P,  and  <r 
but  independent  of  u  and  h,  such  that 

(3.13)  ||u-d.u|l,  _  <  Chlul,  _  V  u  e  lf(n).  0  <  h  <  1. 

n  1 ,  n  L,  n 

Proof.  Since  the  function  u-d,u  is  in  H^CQ),  from  the  Poincare 

n  u 

inequality  we  have 

(3.14)  llu-dj^ull^  <  C(n)  ^  lu-d.j.ul  J  .J.. 

Tee. 

n 

Combining  (3.9),  (3.12),  and  (3.14)  we  get 

4 

(3.15)  ||u-d^u||^^„SC  ^  5  lulf  5  C  /  ;jax  h? 

h 

It  follows  immediately  from  the  definition  of  the  mapping  (2.7)  that 

(3.16)  h;|i  S  max(i,  Dhj.  S  mao<(i,  l)h  V  T  e 

Finally,  estimate  (3.13)  follows  directly  from  (3.15)  amd  (3.16).  ■ 

As  a  consequence  of  the  stability,  approximabi 1 ity,  amd  regularity 
results,  we  obtain  am  estimate  for  the  error  ~ the  H^(n)-norm. 

2 

Theorem  3.4.  For  f  e  L  (11)  let  u  be  the  solution  to  (3.  1)  amd  let  u  be 
the  solution  to  (3.6),  with  S,  defined  in  (3.5).  Then  there  is  a  constamt 
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C  =  C(a,P.<r),  depending  on  a,fi.  and  tr  but  independent 
of  f  and  h,  such  that 

(3.17)  ^  0  <  h  S  1. 

Proof.  It  follows  from  Theorem  3.1  and  standard  results  on  the  approximation 
of  problems  in  variational  form  that 

(3.18)  ||u-u^||^^^<C  inf  llu-,11^ 

Combining  (2.6),  (3.13),  and  (3.18),  and  the  fact  that  u  e  11^(0)  implies 

d,  u  e  S,  ,  we  have 
h  h 

llu-Vi.Q^Chllfllo^^, 

where  C  =  C(a,8,<r).  ■ 

Theorem  3.4  shows  that  the  method  defined  by  (3.6)  is  accurate  auid  robust 
for  the  approximation  of  (3.1),  i.e.,  the  convergence  is  of  first  order  in  the 
mesh  pareimeter  h  with  a  constant  depending  on  a  and  8.  but  otherwise 
independent  of  the  coefficient  a(x).  Thus  the  method  has  the  same  accuracy 
as  the  usual  finite  element  method  based  on  C*^,  piecewise  linear 
approximating  functions  for  smooth  problems. 

Remzirk  3.  1  Approximation  Method  I,  as  we  have  presented  it,  is  baised  on  a 

triangular  mesh.  One  can  also  consider  rectangular  meshes.  Thus  for 

0  <  h  S  1,  let  be  a  peirtition  of  Q  by  rectangles  R  of  diameter  <h 

and  suppose  ^^h^0<hSl  satisfies  a  "minimal  angle  condition"  ((diaun  R/diaun  of 

largest  disk  contained  in  R)  <  <rVR  e  W,  and  V  0  <  h  S  1).  With  !».  we 

h  h 

associate  the  approximating  functions 

=  (V  e  C°(5):  v|^  e  j^|i,  y.  J  VR  € 

V  =  0  on  an}  . 
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The  finite  element  approximation  is  defined  by  (3.6)  with  this 

choice  for  Sj^.  Then  it  is  easily  seen  that  the  airguments  used  to  prove 
Theorem  3.4  yield 

(3.17')  ||u  -  Q  ^  C(a.3,o-)h  llfllg  jj. 

the  S8une  estimate  proved  for  triangular  meshes. 

Remark  3.2.  Method  I  has  an  obvious  one  dimensional  version.  This  one 
dimensional  method  differs  from  the  standard  finite  element  method  based  on 
C^,  piecewise  llneeu'  approximating  functions  in  that  the  coefficient  a(x) 
enters  the  finite  element  calculations  via  its  element -by-element  harmonic 
averages  instead  of  via  its  averages.  It  is  referred  to  as  a  generalized 
displacement  method  (cf.  [4]).  In  the  methods  presented  in  this  paper,  the 
coefficient  a(x,y)  enters  the  calculations  via  various  element -by-element 
harmonic  averages  and  averages,  l.e.,  via  various  element -by-element  moments 
of  l/a(x,y)  and  a(x,y). 

3.2.  Approximation  Method  II. 

In  Method  I  we  chose  shape  functions  that  closely  approximated  the 
unknown  solution.  We  then  used  the  seune  fimctions  for  test  functions,  and  the 
stability  condition  was  Immediate.  In  order  to  ensure  our  methods  were 
conforming,  we  used  curvi linear  triangles.  In  this  subsection,  we  discuss  a 
second  method,  employing  the  trizoigulation  by  ordinary  triangles  shown  in 
Fig.  3.2,  the  trial  functions  used  in  Method  I,  and  piecewise  linear  test 

functions.  Now  the  trial  space  will  be  non-conforming,  but  the  test  space 
will  be  conforming. 

For  h  =  ^,  n  =  2,3,...,  let  be  the  uniform  triangulation  of  Q, 

with  nodes  ~  (ih,Jh),  i,J  =  0 . n,  shown  in  Fig.  3.2. 
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Fig.  3.2 

For  use  in  our  ar  ,]  sis,  we  introduce  the  mesh  dependent  spaces 

(3.19)  H^CQ)  =  {u  €  L^(Q)  :  ul.j.  €  H^(T)  V  T  e 

with  the  norms 


u 


l,h 


u^  dxdy  +  I u I ^  . 

1,  h 

Q 


(3.20) 


u  dxdy  + 


I 


Q 


I  grad  u|  dxdy. 


h 


It  is  clear  that  these  spaces  are  Hilbert  spaces. 

We  define  the  bilinear  form  on  H^(n)xHQ(n)  by 


(3.21) 


U,V)  =  Y  I 
X.-t“  J- 


a  grad  u  •  grad  v  dxdy. 


Tel?.  ’'T 
n 


Clearly  is  bounded  on  H^(n)xHQ(n),  with  a  bovmd  that  is  independent  of 

h.  Moreover,  Bj^(u,v)  =  B(u,v)  V  u,v  e  HQ(n).  Now  we  define  the  trial  space 
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2uid  the  test  space  S, 


(3.22)  S 


l.h  "  { 


V  €  L  (n)  :  V 


€  span{l,J^ 


a(t) 


.y>  V  T  € 


V  is  continuous  at  the  nodes  of 

V  =  0  at  the  boundary  nodes!- 


(3.23) 


S  =  (v  e  C°(n)  :  v|.j.  €  span{l,x.y>.  =  0} . 


We  remark  that  4  Hq  in  general,  so  ^  is  non-conforming  as 

mentioned  above. 

Our  finite  element  approximation  u^^  to  u  is  then  defined  by 

%  ^  ^l.h 

(3.24)  - 


Note  that  the  space  is  not  well  suited  for  a  weak  formulation  of  the 

exact  problem  (3.1).  Nevertheless,  the  error  analysis  of  (3.24)  can  be 
cau’ried  out  in  the  usual  way.  Let  us  suppose  that  a  stability  condition  holds 
for  (3.24),  i.e.,  there  exists  5  =  3(a,^)  such  that 


(3.25) 


inf  sup  |B.  (u,  v)|  >  3(a,B)  >0  V0<h5  1. 

|u|ii,h=i  Ill'll i.n'i 


Since  dim  S,  =  dim  S„  .  ,  (3.25)  implies  that  (3.24)  is  uniquely  solvable. 
1 ,  h  2, h 


For  any  u  e  H^(n^,  we  can  define  P^u  by 


(3.26) 


V  ^  Sl,h 

Bj^(P^u,v)  =  Bj^(u,v)  V  V  e 


It  is  clear  that  P^^  is  a  projection  onto  This  projection  is 

uniformly  bounded  in  h;  in  fact  by  (3.25)  suid  (3.26)  we  have 
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lIPj^ulli  ^  ^  S(a.p)  ^  sup  |Bj^(P^u.v)l  S  C(a,B)llul|j  ^  ^  ^  ^  Hj'^(n)- 

For  u  the  solution  of  (3.1)  and  u^^  the  solution  of  (3.24)  we  have  for  any 
'  ®i,h’ 

"“‘Vi.h  “  Vi.h  “  !(“-):) -Ph(u-*)ll,  „  s  (i*c(.,@)]||u-x||j  1^. 

Thus  we  have  proved  there  exists  a  constant  C  =  C(a,^)  such  that 

(3.27)  llu-Uj^ll^^^.C  inf 

^^^l.h 

(cf.  [1]). 

We  show  now  that  the  stability  condition  (3.25)  holds. 

Theorem  3.5.  There  is  a  positive  constant  6  =  5(a,0)  such  that 

(3.28)  inf  sup  |B  (u,v)|  >  6(a,8)  V  0  <  h  <  1. 


Proof.  Let  a^ 
i .  e. ,  let 

(3.29) 


(0,1) — >1R  denote  the  piecewise  harmonic  average  of  a(x), 


-1 


J 

r^i  1 

dt 

-1 

^  h 

a(t) 

.  J 

X  - 

4 

i-1 


where  =  (x^_^,x^).  For  any  u  €  let  v  e  ^  be  defined  by 


v(P)  =  u(P)  for  all  nodes  P  of  lo.  . 

n 


We  will  verify  now  the  relations: 

(3.30) 


9u  _  dv  du  _  dv 
^dx  ^  dx’  dy  ~  d^  ■ 
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Let  us  first  consider  a  triangle  T  of  the  type  shown  in  Fig.  3.3. 


Fig.  3.3 


On  T  we  have 


+  u(x,  ,y ,) 


x-xi-i 


h 


and 


ulx.yl  = 


r 


(y-yi-i)  _  xi-i 


dt 

a(t ) 


*u(Xj,yj) 


r . 

Xi-i 


r. 

Xl-l 


dt 

a(t) 


dt 
a(t ) 


rs 

Xl-l 


dt 

a(t) 


From  these  two  formulae  we  clearly  get 


^  _  u(xt.y')-u(xi-i.y') 
dx  h 


x,.i 


dt 

a(t) 


dv 

\  ax’ 


dxxdv 

H-  =  H—  in  T 

dy  dy 
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On  the  triangles  of  the  type  shown  In  Fig.  3.4  the  proof  follows  the  same 


X 


Fig.  3.4 


lines.  So  the  relations  (3.30)  aire  proved.  Now  using  (3.30)  and  the  Polncar^ 
Inequality  we  have 

(3.31)  B^Cu.vl  =  ^  f  .  a?H  gjdxdy 

n 

~  =  a|v|^  ^  >  |l|v||^ 

Jf 

n 

To  complete  the  proof  we  still  have  to  bound  ||v||  _  from  below  in  terms  of 

1  » 

llUill^  h"  ^he  relation  (3.30)  we  obtain 


(3.32) 


»  *7’ 


T€^. ’'T 
n 


On  the  kind  of  triangles  T  shown  in  Fig.  3.3  we  have 
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(3.33) 


I u I  dxdy  = 


dt 

r  a(tk  a{t) 

u(x  .y 

i  1  J  I  ^  dt  J  1  J  dt 

j  _  a(t)  J  _  a(t) 


dxdy 


X,-i 


Xi_i 


<  3 


{  ■*■  lu(x^.yj)  l^}dxdy 


T 
2 

=  3  ^  |u(x^,y^)  1^} 


i  C 


I V I ^dxdy 


*  r 


On  triangles  of  the  type  shovm  in  Fig.  3.4  we  have  the  same  estimates. 
Inequalities  (3.32)  and  (3.33)  show  that 


(3.34) 


From  (3.31)  zind  (3.34)  it  follows  immediately  that  (3.28)  holds  with  6  = 

.2 


fcr 

a 

J 


-1/2 


For  u  e  H^(n),  let  d.u  be  the  S,  , -interpolant  of  u,  i.e.,  let 

h  1 ,  h 


dj^u  be  defined  by 


(3.35) 


•V  "  ^l.h 


h 

- 

d,  u 
>•  n 


,  >-(P)  =  u(P)  V  nodes  P  of  lo,  ; 
h  h 


d.u  is  well  defined  since  u  is  well-defined  on  the  nodes  and  since  the 
h 

images  of  the  vertices  of  ainy  T  €  are  nonco linear.  In  the  next  theorem 

we  derive  an  estimate  for  the  interpolation  error  ||u-dj^u||^ 
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Theorem  3.6.  There  is  a  constant  C  =  C(a, ^),  depiending  only  on  a  and  $ 
but  independent  of  u  and  h,  such  that 


(3.36) 


|u-d.u||,  ,  <  Chlul,  _  V  u  e  lf(n),  0  <  h  <  1. 
n  1 ,  n  L»  w 


Proof.  This  proof  is  similar  to  that  of  Theorem  3.3.  Let  u  e  T  € 

aund  Rj.  be  the  smallest  rectangle  containing  T.  Let  T.Rj.  be  the  images  of 
T,R  under  the  mapping  (2.7).  Then,  applying  the  usual  linear  interpolation 
theorem  ais  in  (3.11),  we  have 


(3.37) 


Vl.R,  ^  ‘ 


To  obtain  the  result  in  the  original  variables  we  note  that 


(3.38) 


Ht 


and  (cf.  (3. 10)) 


(3.39) 


'“’V'i.Rt  ^  Vi.Rt- 


By  the  definition  (2.5b)  of  the  semi-norm  1*1,  _  ,  we  have 

L,  Kx 


(3.40) 


''^'2,Rt  “  ''''l,Rt* 


With  (3.38)  -  (3.40),  inequality  (3.37)  yields 


(3.41) 


P 


From  (2.7)  we  have 


(3.42) 


h^^  <  max(  1/a,  1  )hj^^,  >  min(l,l/B)  pj^^. 


Finally,  (3.36)  is  a  consequence  of  (3.41)  auid  (3.42). 


As  a  consequence  of  (3.27),  (3,28),  (3.36),  and  (2.6)  we  obtain 


Theorem  3.7.  For  f  €  L'^(n),  let  u  be  the  solution  to  (3.1)  and  let  u^^ 
be  the  solution  to  (3.24).  Then  there  Is  a  constant  C  =  C(a, 8)  such  that 


(3.43) 


Remark  3.3.  In  Remzu'k  3.1  we  briefly  outlined  Method  I  for  rectangular 
meshes.  Here  we  give  a  rectauigular  mesh  version  of  Method  II.  Let 


h  ’  f  ■  v|„  *  span  {l.  y,  j, 

V  =  0  on  9  n  j- 


V  R  e 


Sg  =  -^v  €  C°(n)  :  vjj^  €  span  |^1,  x,  y,  xy  j-V  R  e  v  =  0  on  9  nj-. 

Our  finite  element  approximation  u^  is  now  defined  by  (3.24)  with  this 
choice  for  and  In  this  situation  we  need  a  hypothesis  on  a(x) 

in  order  to  ensure  stability.  Let 


f  r 

^1-1 

(x)  =  - 

‘  r‘ 


’‘i-i  “■’‘i-i 

♦i'"’  =  IT - T  *1 - — 


and  then  let 


,  =  H-  f  > 


a  dx. 


=  h  ^  <p^  adx. 
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Cj  =  h~^  adx. 


•i-1 


and 


a  dx. 


i-1 


We  assume 


4A^Ci-  (B^.D^r 
4C, 


>  y  >  0,  Vi.h. 


Then  (3.28)  holds  with  5  =  S(a,fi,y)  >  0.  We  therefore  obtain 


(3.43'  ) 


We  remark  that  ^  is  conforming  in  this  rectangular  mesh  case  in  contrast 

with  the  triangular  mesh  caise  in  which  S,  .  is  nonconforming. 

1 » n 


3.3.  Approximation  Method  III. 

In  Method  I  we  introduced  curvilinear  triauigles  in  order  to  ensure  the 
approximating  functions  were  conforming,  while  in  Method  II  we  used  a  special 
triangulation  with  ordinary  triangles  obtaining  a  nonconforming  method.  In 
this  section  we  design  a  conforming  method  based  on  am  arbitrary  triangula¬ 
tion  with  ordinary  triangles. 

For  0  <  h  <  1,  let  1?^  be  a  triauigulation  of  H  by  ordinary  triauigles 
of  diameter  <  h  and  suppose  ^^h^0<h<l  s®^tisfies 

(3.44)  V  Teg.Vh  (minimal  angle  condition) 
and 

(3.45)  ^  £  y  V  T  €  V  h  (quasi-uniform  condition), 

where  h_.p_  have  been  defined  in  (3.3),  (3.4).  Let  P,,...,P  be  the 

11  1  mh 
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nodes  of  loj^.  The  function  denotes  the  usual  piecewise  linear  beisis 

function  associated  with  the  node  Pj  =  J  ~  ^ . i-e., 

is  piecewise  lineau*  with  respect  to  1?^, 


we  remark  that 


I*.  -  >• 


For  each 


J  6  {1 . m^} 


Vj  =  span<0j(x.y),0j(x.y)  0j(x,  y)  (y  -  y^)}  . 


For  the  space  of  approximating  functions  we  choose 


(3.46)  S.  =  {v  : 

h 


•• 

R  :  V  =  ^Vj,  Vj  e  Vj,  v  =  0  on  dfl). 


Our  finite  element  approximation  u^^  to  u  is  now  defined  by 


(3.47) 


%  ^  \ 


[u^.v)  =  L  fv  dxdy  V  V  € 


Uj^  is  the  Ritz  approximation  to  u  determined  by  the  variational  formulation 
(1.3)  and  by  the  space  defined  in  (3.46).  To  study  the  convergence  of 

the  approximation  (3.47),  we  turn  our  attention  to  an  approximation  result  for 

^  Vo<h:Sr 

First  we  show  that  we  can  approximate  u  €  H^(n)  by  a  linear 


combination  of  3d\)’  ^ 


on  Sj,  where  for  J  =  1 . m^^,  Sj  i: 


the  finite  element  star  associated  with  the  node  P  : 

sj 
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S  *  U  T 
J  Telfh 
PjcT 


Let  Rj  be  the  smallest  rectangle  with  sides  parallel  to  the  eoces  containing 
Sj  and  let  Jj  be  three  specific  vertices  of  Rj,  including  all  vertices 
that  lie  on  an.  For  any  u  €  H^(Rj),  we  define  the  span{l,  r  a(t) ’ ^ 


><3 


interpoleuit  of  u  associated  with  Pj  by 


(3.48)  d,  .  u  =  a  +  8 
J.  h 


dt 


f^+ir(y-yj).  (dj  j^u)(P)  =  u(P)  V  P  €  Jj. 


'J 


We  will  prove  the  following  approximability  result. 


Theorem  3.8.  There  is  a  constaint  C  =  C(a, 8)  such  that 

(3.49)  |u-dj^^u|^_S^  .  C(«.8)  lul^  .  Vi  =  0.1.  J  =  1 . m^^. 

u  6  if(n). 

Proof.  Let  j  €  (1 . m^^)  auid  u  e  H^(n)  be  given.  With  the  node  P^  we 

associate  the  finite  element  star  S.  and  the  rectauigle  R..  S.  and  R 

0  %)  0  \) 

are  the  images  under  the  mapping  (2.7)  of  S  and  R  .  Clearly  R.  is  also 

3  3  3 

a  rectangle.  It  follows  from  (2.7)  and  (3.48)  that  d  u  is  the 

J.  h 

span{l.x.y>-interpolant  of  u.  i.e..  d~  u€  spanll.x.y}  and  d.~  u(P)  = 

J  f  J  »  ^ 

u(P)VP€3j.  Thus 


(3.50) 


-  -  ■’'s 

'“"‘^J.h^’i.R  ^  i  ■’•'^*2.8  ’  ^ 

J  J 


Returning  to  the  original  variables  in  (3.50)  (cf.  (3.38)  -  (3.40)),  we  obtain 
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(3.51) 


lu-d  u| 

J  I  O 


C(a.^) 


1  =0,1. 


As  in  (3.42)  we  have 

(3.52)  S  max(i,  2  mln(l,|)pj^^. 

Finally,  (3.49)  is  a  consequence  of  (3.51)  and  (3.52).  ■ 


Before  stating  an  approximation  result  for  ^^h^0<h^l’  prove  a 
technical  result. 


Lemma  3.  1.  Let  (lo.  >_  ...  be  a  family  of  trlansulatlons  satisfvinfi  the 

n  u<nsi 

minimal  angle  condition  (3.44).  Let  P . ,P  denote  the  nodes  of  and 

1  “h  h 

let  S.  be  the  finite  element  star  associated  with  P..  Then  we  caji 

partition  the  set  {P,,...,P  >  of  nodes  into  a  finite  number  of  disjoint 

1  mh 

sets  . I^,  with  I  depending  on  <r  but  independent  of  h,  such  that 

Pj,Pj  €  Ij^,  i  *  J,  implies  n  =  0  (§j  denotes  the  interior  of  S^). 

Proof.  The  proof  is  simple;  in  fact,  we  give  an  algorithm  to  construct  the 

partition.  We  assimilate  the  trlauigulation  to  a  graph,  the  edges  being  arcs. 

Because  of  the  minimal  angle  condition,  a  node  Pj^  hais  a  limited  number  of 

arcs  P,  P,  ,  Q,  =  (P,  :  k  =  l,...,3f,>  being  the  neighbors  of  P, ,  with 

1  1  1|(  1  1 

S  y,  where  y  depends  on  <r  but  is  independent  of  i  and  h.  We  now 

state  the  algorithm.  To  construct  1^  we  do  the  following.  First  take  Pj^ 

in  I,;  then  tedce  the  node  of  smallest  index  s  in  (P . ,P  >\  ({P,>  u 

1  1  mh  1 

Q, ) ,  to  ensure  n  §  =0,  and  so  on  until  the  set  {P,,...,P  }\({P,}  u 

1  Is  1  mh  1 

Q.  >j{P}uQ  u...)  is  empty.  To  construct  I_  we  do  the  same  as  before 

1  S  S  b 

with  the  subgraph  of  nodes  {P . ,P  }\I,.  From  the  minimal  angle  condition 

1  nih  1 

amd  the  construction  of  I^,  a  node  of  this  subgraph  will  have  at  most  y-1 
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sires.  In  this  way  we  construct  .  The  algorithm  will  stop  after  at 

most  T(  steps.  ■ 


Theorem  3.9.  There  is  a  constsuit  C  =  C(a, 0, y.o-),  depending  on  a,p,v,<r 
but  independent  of  u  and  h,  such  that 


(3.53) 


o<hsi. 

J=1 


Proof 


Let  u  €  H^(n)  ajid  let  Ij,...,I^  be  the  pau~tltion  of  the  nodes  of 


given  in  Lemma  3.1.  Then,  since  supp  i/f.  =  S.,  we  have 
n  J  J 

“h  n»h 

(3.54)  l'‘-XVj.h“'?,n= 


‘  'Z  Z 

k=l  j€lk 
I 

-  ‘Z'Z 

k=l  jelk 

t  . 

=  1^  grad(^j(u-dj  dxdy 

k=l  '^n  j€lk 


I 

-z  z  I  (grad  ^j)(u-dj  j^u)  +  grad(u-dj  dxdy 

k=l  Jelk 

IHh  p 

^  {((grad  ^^^J)(u-dJ  j^u)|^+  grad(u-dj  j^u)|^>dxdy. 

J=1  Jsj 

We  note  that  with  the  assumptions  (3.44),  (3.45),  we  have  the  bounds 
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l^jl  i  1  and 


I  grad 


- ; —  £  — I —  < 

min  p_  min  h_  h 

PjcT  PjeT 


Thus  from  (3.54)  we  get 


(3.55)  |u- 


“h  ®h  f 

1=1  1=1  V.  Pj 

j  1  j  1  p  j 


|u-d.  .u|  dxdy 
J  »  *1 


f  2 

+  |grad(u-dj  |  dxdy 
J  ^ 


We  now  use  Theorem  3.8  in  (3.55)  to  get 


(3.56) 


2 

i.a 


mu 

£  2a:(a.8)^ 
J=1 


•4 


min  p_ 
PjsT 


2 

L,Rj 


With  the  assumptions  (3.44)  and  (3.45),  the  following  estimates  are  obvious 


h„  i  2  max  h_  <  2h,  p^  ^  min  p_  >  min  ^ 
PjeT  ^  Pj€T  ^  PjeT 


So  the  inequality  (3.56)  becomes 


(3.57) 


mh  mh 

'“‘Z  Vj.h“'?.Q  - '“'uR,- 

J=1  J=1 


2 

It  remains  to  estimate  V  |u|-  _  .  We  have 

j=i 


(3.58) 


where 


mh 

Z  '“'Ir.  "  Z 


J=1 


Te1?h 


N.J.  =  The  number  of  rectangles  Rj  such  that  T  a  *  0. 
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Let  us  now  show  that  under  the  assumptions  (3.44)  and  (3  45)  the  numbers  N.^. 
caul  be  bounded  independently  of  h  for  all  T  e  Let  T  e  be  given. 

If  T  A  R.  ^  0,  then  P.  lies  within  the  (closed)  disk  D  of  radius 
(v^+l)h  centered  at  the  center  of  T.  To  estimate  the  number  of  nodes  lying 
inside  D,  we  first  estimate  the  number  N  of  triauigles  K  that  lie  inside 
the  disk  D'  of  radius  (v^+2)h.  Since  from  (3.44),  (3.45)  we  get 

k2  w2  2 

^  ^  n — ^  S  Ti  <  area(K), 

4i;  cr  4<r 

we  have  the  estimate 


N 


.22 
4y  O’ 


^  ^  area(K)  S  area(D'  )  =  jr(2+V5)^h^, 


KcD' 


and  hence 


N  i  4i;^<r^(2+V^)^. 


So  for  N.J.  we  have  the  bound 

(3.59)  N_  <  12i;^<r^(2+v^)^  V  T  e  1?^,  0  <  h  <  1. 

T  h 

Finally  combining  (3. 57)-(3. 59)  we  get 

IHh 

<"•“>  '“-Z  Vj.h"'?,!!  ^ 

J=1 

where  C  depends  on  but  not  on  u  nor  on  h.  Inequality  (3.53) 

follows  from  (3.60).  ■ 


As  with  the  Approximation  Method  I,  the  stability  condition  is  immediate 
(cf.  Theorem  3.1).  In  the  saime  way  we  proved  Theorem  3.4,  we  can  prove 

2 

Theorem  3.  10.  For  f  €  L  (£1),  let  u  be  the  solution  to  (3.1)  and  let  u^^ 
be  the  solution  to  (3.47).  Then  there  is  a  constant  C  =  C(a,R,y,(r)  such 
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that 


(3.61) 


Note  that  in  the  proof  of  Theorem  3. 10  we  use  the  fact  that  u  e  H^dl) 
mh 

implies  ^  h'^  ^  This  is  true  because  Jj  contains  any  vertices 

J=1 


of  Rj  that  lie  on  50. 

3.4  Comments  on  Methods  I.  II.  Ill 

We  have  described  three  methods  for  approximating  the  solutions  of 

problems  of  the  type  depicted  in  Fig.  1.1.  The  usual  finite  element  method 

1+c 

is  inaccurate  for  these  problems  since  the  solutions  may  not  be  in  H  (0) 
for  any  e  >  0. 

Methods  I  and  II  are  closely  related.  The  central  idea  in  these  methods 
is  to  exploit  the  existence  of  a  mapping  from  the  general  element  to  the 
reference  element  that  transforms  the  special  shape  functions  into  polynomials 
auid  the  unknown  solution  into  a  smooth  function,  and  thereby  obtain  a  good 
convergence  rate.  For  singular  corner  behavior  aind  homogeneous  material,  this 
idea  is  exploited  in  [5]. 

It  is  advantageous  to  use  rectangular  meshes  in  Q  that  Eire  aligned  with 
the  direction  of  the  unidirectional  composite,  as  described  by  a(x,y), 
because  they  are  the  images  of  rectEuigular  meshes  on  0.  The  major  difference 
between  Methods  I  and  II  is  in  their  treatment  of  the  right  hEUid  side  f. 

Since  with  Method  II,  f  enters  the  computation  through  integrals  of  f  times 
the  usual  piecewise  linear  test  functions  (eis  opposed  to  integrals  of  f 
times  the  special  test  functions  (cf.  (3.5)),  Method  II  is  preferable  when 
many  right  hsuid  sides  must  be  treated.  On  the  other  hand.  Method  II  is  less 
stable  tham  Method  I,  leading  to  larger  constants  in  the  error  estimates  (cf. 
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(3.17),  (3.17'),  (3.43),  (3.43')).  In  fact,  for  quadrilateral  meshes  Method 
II  may  not  converge  for  some  a*s;  see  the  hypothesis  on  a(x)  in  Remark  3.3. 
We  note  that  for  trieuigular  meshes  Method  II  always  converges. 

Method  III,  although  similar  in  its  use  of  good  local  approximating 
functions,  (e.g. ,  functions  satisfying  the  differential  equation)  has  a 
rather  different  character  than  Methods  I  and  II.  In  Method  III  the  alignment 
of  the  mesh  does  not  play  a  role.  Finite  element  approximating  spaces  based 
on  shape  functions  satisfying  the  differential  equation  have  been  suggested 
and  employed  in  various  contexts.  The  main  problem  in  their  use  is  the 
enforcement  of  some  type  of  conformity.  This  can  be  done  by  various  hybrid 
methods,  e.g.  (see  [2],  [12],  e.g.).  These  are,  however,  problems  in  ensuring 
the  stability  of  these  methods,  and  some  of  these  problems  have  not  been 
satisfactorily  resolved.  In  contrast.  Method  III  proposed  here  has  no 
problems  of  this  type,  and  is  very  accurate  and  robust.  For  some 
computational  aspects  of  a  similar  method  employing  harmonic  polynomials  in  a 
p-version  fashion  and  applied  to  the  solution  of  Laplace’s  equation  see  [15]. 
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4.  Methods  for  Problems  with  Curvilinear  Unidirectional  Coefficients 


The  methods  presented  in  Section  3  cover  problems  on  rcct^ngulc^'  domains 
with  coefficients  that  globally  vary  sharply  in  at  most  one  direction,  i.e., 
that  are  straight  line  unidirectional.  Here  we  extend  Methods  I  and  III  to 
cover  coefficients  that  locally  vary  sharply  in  at  most  one  direction,  i.e., 
that  are  curvilinear  unidirectional,  and  to  cover  domains  with  curved 
boundEU'ies.  Method  I',  the  extension  of  Method  I,  will  be  based  on  quadri¬ 
lateral  and  triangular  elements  and  Method  III',  the  extension  of  Method  III, 
will  be  beised  on  triangular  elements. 

4.1.  Method  I' 

Consider  the  boundary  value  problem  (1.1)  and  suppose 

•  for  1  5  i  <  n' ,  )  is  am  open  subset  of  n  and  a 

coordinate  system  satisfying  conditions  (i)-(vi)  in  Section  2  and  for 

n'  +  1  5  i  ^  n,  where  n'  <  n,  is  an  open  subset  of  n 

and  a  coordinate  system  satisfying  conditions  (i)-(iii),  (v),  (vi)  in 

Section  2,  i.e.,  C  =  =  ?.(x,y),  t)  =  t).  =  Ti.(x,y)  and  if  (x,y) 

12  1 

ranges  over  n, ,  then  (?, ,  tj  )  ranges  over  H, '  =  (^_  ,?„  )  x  (7)„  , 

III  I  \i\ 

2 

7J-),  where  ,  ^,(x,y),  ii.(x,y)  satisfy  conditions  (i)-(iv),  (vi) 

Wi  II  I 

if  i  <  n'  amd  conditions  (i)-(iii),  (vi)  if  i  >  n'  +  1  (let  E. 

1 

denote  the  union  of  the  interior  edges  of  ) ; 

•  covers  n  in  the  sense  that 

(4.  1)  n  = 

and 

(4.2)  an  =  Interior  (n^  n  an)  In  an); 

•  for  1  $  i  5  n' ,  we  have 

(4.3)  a(x,y)  =  a(x^(^,7j),  yj(C,T)))  =  a^(^)  V  (x,y)  €  n^, 
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where  x  =  y  =  is  the  inverse  of  the  mapping 

?^(x,y),  T}j^(x,y)  (a^  here  plays  the  role  of  a'  in  (2.12)),  and 

for  n'  +  1  <  i  <  n,  a(x,y)  is  smooth  on 

With  {  satisfying  these  aissumptions,  for  each  i  let  0^ 

be  the  result  of  pulling  each  interior  edge  of  Ilj^  (cf.  condition  (vi)  of 

Section  2)  a  distance  d  towards  the  center  of  Then  the  are  open 

sets  of  the  type  considered  in  Theorem  2.3  (i.e.,  0.  c  £1.,  0.  cc  £1.  if 

1111 

£2.  n  3£2  =  0  and  dO.  n  a£2.  c  3£2  if  Q.  n  3£2  ^  0)  and  {0.}?  ,  satisfies 
1  111  1  1=1 

(4.1)-(4.2).  provided  d  is  sufficiently  small.  Note  that  d  =  dist{0^,E^). 
We  consider  d  to  be  fixed. 

We  note  that  if  (1.1)  corresponds  to  problems  of  the  type  depicted  in 
Figs.  1. 1-1.3  or  to  a  smooth  interface  problem  modeled  as  a  composite 
material,  then  {(£1^,  can  be  chosen  to  satisfy  the  conditions 

outlined  above. 

We  note  that  these  aissumptions  imply  that  3  £1  is  a  piecewise  smooth 

2 

(C  )  curve  with  vertices  with  suigular  measure  a  satisfying  0  <  a  <  tt;  in 
particular  f.  has  no  reentrant  vertices. 

Remark  4. 1.  If  our  proclem  satisfies  these  assumptions  we  say  that  a(x,y) 
locally  varies  sharply  in  at  most  one  direction.  Such  coefficients  aire,  as 
indicated  earlier,  also  called  (curvilinear)  unidirectional. 

With  n,  £2j . £2^,  n' ,  and  a(x,y)  satisfying  the  hypothesis 

described  above,  we  now  describe  the  meshes  we  will  employ.  For  0  <  h  <  1, 
let  =  (T)  be  a  mesh  on  £2  consisting  of  curvilinear  (closed)  quadri¬ 
laterals  or  triangles,  and  satisfying  he  following  properties: 

•  Each  T  is  contained  in  some  0^  :  T  c  i  -  i(T)  <  n; 

•  If  1(T)  <  n' ,  then  T  is  the  image  of  a  rectangle  T'  in 


under  the  mapping  x  =  y  =  (-j-)  (C.  t^)  . 


i.  e.  , 


T  =  {(x,y)  ;  X  =  x^^^j(C.i?),  y  = 


i  ^  ^  ,7)1  5  T)1  i  TJ  S  7|^  <  7)^  }. 

Ukt)  t  t  wkt)  “i(t)  t  t  Hkt) 


where 


(4. 4a) 


(4.4b) 


IC^  -  ^  h,  IT,^  -  7,^1  <  h. 


1€^  - 

-1  ^ 


where  1  <  <r  <  «  is  independent  of  the  mesh.  The  mapping 

12  12 

7j^.j.j)  maps  T  onto  T'  ®  (?.j.,  ?.j.)  x  (tj.^.,  tj.^.)  and  T'  is  mapped 
onto  the  reference  rectangle  =  (0,1)  x  (0,1)  by  the  mapping 


(4.5)  t  =  I, 


.1  ^i(T)  J-i  ^i(T)  J  1  ^i(T) 

T _  _  ^KT) _ (T) 


.  7)  =  71 


.1  ^i(T) 


•Cj 

c 

.^1  ^il 


r}  -  Vj 

i(T)  ^  ~2  T 

Vj  -  Vj 


Thus  the  composition  of  these  two  mappings  maps  T  onto  T*^  =  T°, 
and  the  inverse,  F.^.,  of  the  composition  maps  T^  onto  T. 

If  i(T)  2  n'  +1,  then  T  is  the  image  of 


T^  =  j^T^,  if  T  is  a  quadrilateral 


T^  =  a  reference  triangle,  if  T  is  a  trieingle 
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(4.6) 
8Lnd 

(4.7) 1J 


under  a  mapping  F.^.  satisfying 


F.p  is  invertible,  and  F.^.  euid  F.j.  are  smooth, 
I^tK.oo.T®  “  ^  K,co,T  ^  ’ 


Ft'o.w.T® 


s  suplJ  (?.?i)l  5  Ch^,  IJ  .1 

r?  rj 

(?.?i)€TO 


O.oo.T 


a  sup  IJ  _i(x,y)  1 
(x,y)er 


5  Ch’2, 


where 


|G| 


I,  00,  Q 


=  sup  l|D^G(t,s)||  ,  2  n2i 

(t,s)€Q  ’ 


|D^G(t,s)||  2  02^  =  max  ||D^G(t,  s)  (y, . 

^  y.€R2  1 


yieR2 

II  n  11^1 

l<i<^ 


•  . 7^)  I 


and 


2 

11*11  =  the  Euclidean  vector  norm  on  R  , 


J|,(t,s)  =  Jacobian  of  G  at  (t,s). 


The  constant  C  in  these  estimates  is  independent  of  the  mesh.  We 
easily  see  that  the  mapping  F.^.  :  T°— >T  defined  above  for  i(T)  < 
n'  satisfies  parallel  aissumptions.  Hence  we  have  T  =  F.j.(T^)  for 
all  T,  and  it  is  convenient  to  associate  the  mesh  ^  =  {T}  with  the 
set  of  mappings  {F.^.}. 

The  stauidard  compatibility  condition  is  satisfied.  Suppose  that  T^ 
and  T_  are  quadrilaterals  with  a  common  edge  t  :  t  =  T,  nT_. 

See  Fig.  4. 1;  note  that  we  are  using  two  copies  of  the  reference 
recteuigle  T^. 
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Fig.  4.1 

Assuming  that  t  Is  the  Image  of  the  vertical  line  segment  {(tiTj)  : 

1=1, 01^51}  under  both  F_  and  F-  ,  we  require  that 

Tl  12 

(4.8)  F-  (1,7?)  =  F._  (I.tJ)  .0  <  tj  <  1. 

ll  12 

If  I  Is  the  Image  under  F_  of  a  different  edge  of  the  reference 

12 

rectangle,  we  would  modify  (4.8)  In  aji  obvious  mauiner.  Also,  if 
either  or  is  a  triangle,  the  compatibility  condition  would 

be  modified  in  an  obvious  way. 

We  point  out  that  our  mesh  matches  the  (curved)  boundary  of  (2  by  meeins 
of  blending  (non- iso parametric)  elements. 

Remark  4.2.  In  the  quadrilateral  element  case,  verification  of  (4.6)  sind 
(4.7)  usually  proceeds  along  the  following  lines.  Let  T*  denote  the 
straight  line  quadrilateral  with  vertices  “  1,2, 3, 4,  coinciding  with 
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those  of  T  (see  Fig.  4.2),  let 
iij.  =  diam  T*. 

=  diajn  of  largest  disk  contained  in  T*, 

=  max  {  |cos  { "  ^i  ^  *  ^^i-1  l^iS4  (mod  4)} 


Fig.  4.2 


One  then  medces  assumptions  on  the  perturbation  ♦  that  imply  (4.6)  and  (4.7) 
are  satisfied.  This  procedure  is  outlined  for  isoparametric  quadrilateral 
elements  in  [7,  Exercise  4.3.9].  The  parallel  procedure  for  triangular 
isoparametric  finite  elements  is  carried  out  in  [7,  Theorem  4.3.3]. 
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It  remains  to  describe  our  shape  functions.  On  T  e  we  use  the  shape 
fxinctions 


(4.9) 


1.  F"^j(x.y),  F‘^2(x.y).  F"^^  (x.y) 


a  quadrilateral 


and 


(4.9')  1.  F.j.^j  (x.y),  F.j.^2(^*y)»  T  is  a  triangle, 

where  F.^.^  (x,y)  =  (F.^.^^  (x.y) ,  F.j.^2  y) )  •  i.e.,  we  use  the  pull-back  poly¬ 
nomials  determined  by  the  bilinear  shape  functions  1,  i),  |t)  in  the  quad¬ 
rilateral  case  and  by  the  line2Lr  shape  functions  1,  t)  in  the  trieingulau’ 

case.  For  i(T)  i  n'  we  easily  see  that  the  functions  in  (4.9)  su'e 


(4.10)  1, 


pi  ^i(T) 


pi  ^i(T) 


(x,y)-i).^ 


1  J 


-  T), 


.€i(T)(x,y) 


ri 


dt 

^i(T) 


dt 

^i(T) 


Then  we  let 
(4. 11) 


S,.  =  {v  € 


L^(n) 


e  spaui  of  the  shape  functions  on  T, 


V  is  continuous  at  the  nodes  of 

h 

V  =  0  at  the  boundary  nodes) . 


Because  of  the  above  assumptions,  in  pau'ticuleu'  (4.8),  we  see  that  c 
11^(0),  i.e.,  is  conforming.  The  Interpol  ant  of  u  is  defined  by 


(4. 12) 


’'‘h  ^  ^h 

(d  u)(P)  =  u(P)  V  nodes  P  of  . 
'■  h  h 
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Because  of  our  choice  of  shape  functions,  dj^u  is  a  good  approximation  to  u. 

In  Fig.  4.3  we  show  a  typical  part  of  the  mesh  on  fl.  '^e  show  the  sets 
and  as  well  as  the  elements  of  the  mesh.  Note  thai  ,  bh  the 

auid  the  elements  fit  the  geometry  of  the  fibers.  In  Fig.  4.4  we  show  the 
mesh  in  a  neighborhood  of  the  boundary  of  fl.  We  see  in  particular  the 
interior  and  the  boundary  edges  of  the  In  Fig.  4.5  we  show  a  typical 

mesh.  We  do  not  show  the  sets  and  0^,  but  do  show  the  areas  where  the 

coefficient  a(x,y)  is  smooth  and  where  it  is  rough  (the  areas  with  the 
fibers).  Note  that  in  the  area  of  the  fibers  we  use  quadrilaterals  elements 
while  in  the  area  where  a(x,y)  is  smooth  we  use  both  quadrilateral  aind 
triangular  elements.  Obviously  triangular  elements  cannot  be  avoided,  but 
quadrilateral  elements  are  preferable  because  they  usually  lead  to  higher 
accuracy  (although  with  the  sane  rate  of  convergence). 


-  Boundary  of  fij 

- Boundary  of  O, 

.  Boundary  of  the  Elements 


Fig.  4.3.  Typical  configuration  of  the  sets  £1^,  0^, 
and  the  elements  inside  £1. 
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The  approximation  property  of  the  spaces  is  formalized  in 


Theorem  4.1.  There  is  a  constant  C  =  C(a, /3,  tj^  ) . ^^n'’^n'^’ 

(n  .  a|  .  a|_  ,  d)  depending  on  a.  . 

n  +1  n  In^*  1  1 

^n'+l’ ■  ■  ■ ’^n’  ^|n  ^|q  ’  independent  of  h  and 

n'  +  1  n 

u,  such  that 


(4. 13) 


Proof.  Consider  T  €  and  let  dj.u  be  defined  by 

d.j,u  €  span  {shape  functions  on  T> , 

(d.j.u)(P)  =  u(P),  for  all  vertices  P  of  T 
(cf(4.12)).  For  i(T)  S  n'  we  see  that  d.j.u  is  well  defined  by  noting 
successively  that  (d.pU)'  (where  the  prime  denotes  the  transformation  from 
the  variables  (x.y)  to  the  variables 


span  •< 


•^i(T) 

dt 

,  -1  a. 1  1 

?T  i(T)  -  TJ.^  Vj 

-2  ’  2  1  ’  2  1 
Aj  Vj 

.pi  ^i(T) 

?T 


A(T) 

dt 


^i(T) 


pi  ^i(T) 


-interpolant 


of  u'  at  the  points  P' ,  that  (d.j,u)'  (where  the  tilde  denotes  the  trans¬ 
formation  from  the  variables  (j)  *  (X)  ^  variables 

r^i(T)  dt 


(4. 14) 


?i 


(T) 


i(T) 


'i(T) 


=  7,) 


’Hi 


is  the  spand,  ?^x)’  ^l(T)’ 


at  the  points 
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P' ,  that  u'  €  from  Theorem  2.3,  emd  that  the  points  P' ,  for  P  a 

vertex  of  T,  form  the  vertices  of  a  rectamgle  in  Note  that  the 

variables  ?jcp)  ^i(T)  here  been  defined  differently  than  in 

(4.5).  This  is  necessary  in  order  that  the  set  depends  only  on  i(T) 

aind  not  on  T.  For  i(T)  i  n'  +  1  we  see  that  d^u  is  well-defined  by 
noting  successively  that  d~u  (where  the  tilde  denotes  the  transformation 
from  the  variables  (x.y)  to  the  variables  ?  =  Ft  ,(x,y),  t)  =  F._  _(x,y))  is 
the  span{l,  7),  |Ti}-interpolant  of  u  if  T  is  a  quadrilateral  and  the 
span{  1,  7?} "Interpol ant  of  u  if  T  is  a  triangle,  that  u  6  H^(T*^),  from 

standard  elliptic  regularity  results  since  T  c  smd  a(x,y)  is  smooth 

on  and  that  the  pxaints  P  aire  the  vertices  of  the  rectangle  T*^.  We 

note  that  these  observations  show  that  d.u  in  (4.12)  is  well-defined. 

n 

Now  we  estimate  |u  -  d.j.u|j^  j.  First  suppose  i(T)  <  n' .  Changing 
variables  we  obtain 


|u-  d.j.u| 


2 

l.T 


/•  f 

'  a(u  -  d.^u) 

2 

a(u  -  d.j.u) 

Jf  ■ 

ax 

ay 

] 

dx  dy 


a(u'  -  (d.j,u)') 

J  •J’  /  ^ 

^^i(T) 

«UT) 

a(u'  -  (d_u)') 

1  I^' 

+ 

X 

^’’i(T) 

Is'-Si  ’>1(1)1 

a(x,  y) 


■  -7 - - - 

^[^i(T).^i(T) 


<  C 


■ 

a(u'  -  (d.j.u)') 

2 

a(u'  -  (dLj.u)') 

2, 

•^T' 

^^i(T) 

^"KT) 
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^  lacG'  -  (d~u)')|^ 


13(11'  -  {d~u)')|^ 


If'  ^^i(T) 


'^^i(T)‘^^i(T) 


^  C  (“•P*?i(X)*’>l(T)^  |u'  dju'l 

where  d^u'  =  (d^u)'  denotes  the  span{l,  ^^(x)'  ^i(T)^i(T)^” 

interpolant  of  u'  on  T' .  Here  ?^x)’  ^i(T)  ^  defined  in  (4.14). 

Thus  by  stamdard  approximation  results  for  bilinear  functions  (cf.  Theorem 


3.1.4  in  [7])  we  have 


l“- Vl.T  ^ 


and  hence,  for  1  <  J  ^  n' , 


(4. 15a) 


Z  Vl.T  ‘ 


i(T)=J 


I  '“''b' 


i(T)=J 


SCh2|u'|2^g^ 


where  C  =  C(a, 0,  (^^, ) . ^^n'’^n'^^’  consider  i(T)  >  n'  +  1.  Using 

(4.6)  auid  (4.7)  and  the  usual  proof  of  approximation  results  (cf.  proof  of 
Theorem  4.3.4  in  [7]),  we  obtain 


Vl.T  ^  ^ 


zind  hence,  for  n'  +  1  S  j  ^  n. 


(4. 15b) 


Z  lu-VI^  T^Ch"  ^  llull"  ,^  S  Ch^llull"^ 
T1=1  i(T)=J 


i(T)=J 
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For  any 


T,  (d,u) 
n 


=  d^u,  and  thus  from  (4.15a,b)  we  have 


(4.16) 

-Z  Z'“'V?,T*  Z  Z  '“'vIt 

l<J<n'  i(T)=j  n'+l<J<n  i(T)=j 


^“"(Z  '“'to  ,*  Z  '“"to,]' 

l<J<n'  n'+l<J<n 

where  C  *  C(a,^,  (?,,  ,,7)  ,).  From  Theorem  2.3  we  have 

11  n  n 

(4. 17)  X  • 

l<J<n' 


where  C  =  C(a,0,  . a(x,y)  is  smooth  on 

Jlj  for  j  >  n'  +1,  from  standard  elliptic  regularity  results  we  have 

(4.18)  ^  '<0.  ^W<Q. 

n'+15J<n 

where  C  =  , . . . ^Ijj  . ^  Ijj  »  As  a  direct  consequence  of 

n'  +  l  n 

(4. 16)-(4. 18)  we  get  (4.13),  as  desired.  ■ 


Our  finite  element  approximation  u^  to  u  is  now  defined  by 


(4. 19) 


% 

B( 


€ 


V''- » 4 


fvdx  V  V  €  S.  . 

h 


Since  we  are  using  Sj^  for  both  the  test  and  trial  space,  stability  is 
Immediate.  Approximabi lity  hsis  been  established  in  Theorem  4.1.  We  thus 
have 
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Theorem  4.2.  Suppose  £2,  . ^  satisfy  the 

assumptions  in  the  first  part  of  this  subsection.  Suppose  u  is  the  solution 


of  (1.1)  and 


is  the  solution  of  (4.19).  Then  there  is  a  constant  C  = 


C(a,3.(r.  (?  Tj^) . V+1 . ”n’  ^|n  ,  , . ^|j2  ’ 

n  +1  n 


d)  such  that 


(4.20) 


I''  “  %h,n  ^  ^  ^  ^  0  <  h  <  hg. 


4.2.  Method  II T 

Consider  the  boundary  value  problem  (1.1)  and  suppose  ((£2^, 

{(£2.,  ^ +  1  •  ^  described  in  Subsection 

4.1.  For  0  <  h  <  h^  let  lo^  =  {T}  be  a  triangulation  of  Q  by  ordinary 

triangles  together  with  curvilinear  triangles  which  fit  the  curved  part  of 

3£2,  all  of  diameter  <  h.  For  any  T  €  lo,  let  T*  be  the  ordinary  triangle 

h 

with  the  same  vertices  as  T.  Then  =  {T*}  is  a  triamgulation  of  £2  by 

ordinary  triangles,  but  U  T*  =  U  T*  is  a  polygonal  approximation  to  £2  and 

T‘6i?*  T€)?h 

h  " 

not  an  exact  fit  of  £2.  We  assume  all  T*  have  diameter  <  h  and  that 

...  satisfies  the  minimal  angle  condition  (3.44)  and  the  quasi-uniform 
h  o<hSho 

condition  (3.45).  Let  (P.  =  (x ., y  ,)>'"*',  be  the  noJes  of  and  let  \(i . 

J  J  J  J=1  h 

denote  the  piecewise  linear  basis  function  corresponding  to  P.  (and  the 

triangulation  T  is  curvilinean,  then  by  restricting  the  domain  of 

definition  of  or  by  extending  i/i  an  a  linear  function  we  can  assume  ip . 

3  3  3 

is  linear  on  T,  and  hence  that  ^  is  continuous  on  £2  and  linear  on  each 

T.  Let  S  =  U  T  be  the  finite  element  stan  ansociated  with  P..  Now  it  is 
■J  T€l?h 


easily  seen  that  if  h^  is  sufficiently  small,  then  any  will  lie  in  some 

0^:  Sj  c  1  ^  i(J)  ^  n.  Let  c  0'^^^  be  the  image  of  n  £2 

i  (  1 ) ' 

under  the  mapping  ^  smallest  rectangle  with 

sides  pairallel  to  the  axes  containing  ,  and  let  be  the  preimage 

\)  \J 
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of  rJ'J’'  under  Ruj,).  since  Isa 

rectauigle,  and  hence  lies  in  Define  ^  three 

specific  vertices  of  including  all  vertices  that  lie  on  5  11. 

Next  we  define  our  space  of  approximating  functions.  For  J  =  l,...,m^ 


(4.21) 

Vj  =  =  span 


.Cl  (x.y) 


Lj(x.y),  #j(x.y)  ^jla.y)  (njU.y)  -  ,J_l| 


if  1  <  i  =  i(j)  <  n' 


i|^^j(x,y),  i/»j(x,y)  (x-Xj),  i^j(x.y)  (y-yj)| 


if  n'  +  1  <  i  =  i(j)  <  n  3md  tl^  f|  an  =  z 


span|i/»j(x,y),  i/(j(x,y)  ICj (x.  y)-?^^  ] ,  ^j(x,y)  [ (x,  y)-T}^j| 


if  n'  +  1  ^  i  =  i(J)  <  n  and  n.  p) 


an  *  0. 


The  3rd  line  in  this  definition  has  been  stated  for  the  case  in  which  the 

11  JZ  I 

preimages  of  the  points  (^_  ,  7j„  ),  (CT  7?_  )  lie  on  an.  In  other 

Wj ,  Wt 

situations  we  could  modify  the  definition  in  an  obvious  way.  Then  for  the 
space  of  approximating  functions  we  choose 


(4.22) 


S.  =  (v  :  n-^R 
h 


’ " '  Z^j’  ''j '  ''j' 


V  =  0  on  an) . 


Our  finite  element  approximation  u^  to  u  is  now  defined  by 


(4.23) 


lUj^.v)  =  J  fv  dxdy  V  v  €  Sj^. 
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Since  we  are  using  the  space  for  both  the  test  and  the  triai  space, 

stability  is  immediate.  To  study  the  convergence  of  the  method  (4.23)  we  need 
to  prove  am  approximation  result  for  the  spaces  ^^h^O<h<l'  done  by 

combining  the  ideas  of  Subsections  3.3  and  4.1. 

Let  d^^jl^u  be  defined  by 


(4.24) 


span 


(x.y) 


1. 


dt 

aTTET 


.  T)^(x,y) 


if  1  <  i  =  i(J)  <  n' 


spanj^l,  x-Xj,  y'yjj"  If  n'+l<i  =  i(j)<n 


and  n.  n  an  =  0 
1 


span|l,  ?^(x,y)  -  y)  - 


if  n'  +  15i  =  i(J)<n  and  £1^^  n  an  =  0. 


(d^K^u)  (P)  =  u(P)  V  P  e 
^  J.h  j 

dj^jl^u  is  a  good  approximation  to  u  on  S^,  as  made  precise  in 


Theorem  4.3. 


There  is  a  constant  C  =  C(a,  0,  (C, .  tj,  ) . (f  ,7)  ))  such  that 

11  n  n 


(4.25)  |u 


d^/;J’u|  ^ 

J,  h  k,  Sj 


h^UJ) 


c  -  lul  i(J)  .  if  1  <  i(J)  <  n' 

P^(J)  ’  J 


‘i 


h!l(J) 

C  -gJ -  Hull  Ri(J)  .  if  n'+l  <  i(j)  <  n  , 

P^(J)  ’  J 

Rj 


for  J  =  1 . mj^  .  k  =  0.1. 
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We  omit  the  proof  of  this  result  since  it  is  similar  to  that  of  Theorem  3.8. 


The  approximabi  1  ity  result  for  {S. 

n  0<nSh, 


is  given  in 


Theorem  4.4. 


that 

(4.26) 


There  is  a  constamt  C  =  C(a, 0,w,<r,  . 


n>h 

.  n'  n 

"  Z  "“"2.0, 

k=l  i=n'+l  ■> 

J=1 

such 


We  omit  the  proof  of  this  result  since  it  is  similar  to  that  of  Theorem  3.9. 

Finally  ais  a  consequence  of  (4.26),  Theorem  2.3,  auid  standau'd  elliptic 
regularity  results  we  have 


Theorem  4.5.  Suppose  £1,  . (€  ,t)  ),  n' ,  auid  a(x,y)  satisfy  the 

11  n  n 

assumptions  in  the  first  part  of  this  subsection.  Suppose  u  is  the  solution 
of  (1.1)  aind  u^  is  the  solution  of  (4.23).  Then  there  is  a  constant  C  = 
C(a,6,y,{r,  (^^,7)^), .  .  .  ,  (^^,7j^),  n' ,  d)  such  that 

(4.27)  n  -  ChilfllQ  Q  V  f  e  L^(n),  0  <  h  <  h^. 


4.3.  Comments  on  Method  I'  auid  III' 

The  diffe”ences  and  similarities  of  Methods  I'  auid  III'  are  similair  to 
those  of  Methods  I  and  III,  which  were  discussed  in  Subsection  3.4.  We  note 
that  with  Method  I'  we  have  to  fit  the  elements  to  the  geometry  of  the  fibers 
of  the  composite,  ais  seen  in  Fig.  4.5.  This  is  not  necessary  in  the  case  of 
Method  III',  and  this  freedom  could  be  utilized  in  many  situations.  For 
exaunple,  suppose  the  coefficient  is  changing  rapidly  but  not  abruptly  along  a 
line.  Then  Method  III'  could  be  used,  leading  to  an  enrichment  of  th“  usual 
finite  element  space  by  special  shape  functions  in  the  neighborhood  of  the 
line. 
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Implementational  considerations  and  computational  studies  of  Methods  I, 


II,  III,  I',  and  III'  will  be  presented  in  a  forthcoming  paper. 


References 

1.  Babu§ka,  I.,  Error  bounds  for  finite  element  method,  Numer.  Math.  16 
(1971),  pp.  322-333. 

2.  BabuSka,  I.,  Oden,  J.T.,  and  Lee,  J.K. ,  Mixed-hybrid  element  approxima¬ 
tion  of  second-order  elliptic  boundary  value  problems.  Part  2,  Weak 
hybrid  methods.  Comp.  Meth.  Appl.  Mech.  Engrg.  14  (1978),  pp.  1-22. 

3.  BabuSka,  I.  and  Osborn  J. ,  Numerical  treatment  of  eigenvalue  problems  for 
differential  equations  with  discontinuous  coefficients.  Math.  Comp.  32 
1978,  pp.  991-1023. 

4.  BabuSka,  I.  and  Osborn  J. ,  Generalized  finite  element  methods:  Their 
performance  and  their  relation  to  mixed  methods,  SIAM  J.  Numer.  Anal.  20 
(1983),  pp  510-536. 

5.  BabuSka,  I.  and  Osborn  J. ,  Finite  element  methods  for  the  solution  of 
problems  with  rough  data,  in  Singularities  and  Constructive  Methods  for 
Their  Treatment ,  Springer-Verlag  Lecture  Notes  in  Mathematics  1121, 
Springer-Verlag,  Berlin-New  York,  1985  (P.  Grisvard,  W.  Wendland,  and 
J.R.  Whiteman,  Editors),  pp.  1-18. 

6.  Bernstein,  S. N. ,  Sur  la  generalisation  du  probieme  de  Dirichlet,  Math. 
Ann.  62  (1906),  pp.  253-271,  69  (1910),  pp.  82-136. 

7.  Ciarlet,  P. ,  The  Finite  Element  Methods  for  Elliptic  Problems,  North 
Holland,  New  York,  1978. 

8.  Ciaurlet,  P.G.  ,  Natterer,  F.  ,  and  Vau'ga,  R.  S.  ,  Numerical  methods  of 
higher-order  accuracy  for  singular  nonlinear  boundary  value  problems, 
Numer.  Math.  15  (1971),  pp.  87-89. 

9.  Crouzeix,  M.  and  Thomas,  J.M. ,  Elements  finis  et  probiemes  elliptiques 
degeneres,  rairo  r-3,  i  (1973),  pp.  77-104. 

10.  Gordon,  W. J.  and  Hall,  Ch.  A.,  Transfinite  element  methods:  blending 
function  interpolation  over  arbitrary  curved  element  domains,  Numer. 
Math.,  21  (1973),  pp.  109-129. 

11.  Grisvard,  P. ,  Elliptic  Problems  in  Non-smooth  Domains,  Pitman  Advanced 
Publishing  Program,  Boston-London-Mel bourne,  1985. 

12.  JlrouSek,  J. ,  Venkatesh,  A. ,  Hybrid-Treff tz  plane  elasticity  elements 
with  p-method  capabilities,  Int.  J.  Num.  Methods  Engrg.  (to  appear). 


67 


13.  Ladyzhenskaya,  0.  and  Ural’tseva,  N.  ,  Linear  and  Quasilinear  Elliptic 
Equations,  Academic  Press,  New  York,  1968. 

14.  Lions,  J.L.  aind  Magenes,  E.  ,  Non- homogeneous  Boundary  Value  Problems  and 
Applications,  Vol.  I,  Springer-Verlag,  New  York-Heidelberg-Berl in,  1972. 

15.  Melenk,  J.M.  ,  Finite  element  methods  with  harmonic  shap>e  functions  for 
solving  Laplace’s  equation  (to  appear). 

16.  Moussaoui,  M.  auid  Ziani,  A..  Sur  1’ approximation  des  solutions  de 
certains  types  d’6quations  diff6rentielles  ou  aux  d6riv6es  partielles  A 
coefficients  non  rAguliers,  Annales  de  I’Enit,  4  (1990),  pp.  9-35. 

17.  Osborn,  J. ,  The  numerical  solution  of  differential  equations  with  rough 
coefficients,  in  Advances  in  Computer  Methods  for  Partial  Differential 
Equations  -  IV,  Proceedings  of  the  Fourth  IMACS  International  Symposium 
on  Computer  Methods  for  Partial  Differential  Equations  held  at  Lehigh 
University,  1981  (R.  Vichnevetsky  and  R.  Steplemaii,  Editors),  pp.  9-13. 

18.  Osborn,  J. ,  Generalized  finite  element  methods,  in  Proceedings  of  10th 
IMACS  World  Congress  on  Systems  Simulations  and  Scientific  Computation, 
held  in  Montreal,  Canada,  1982  (R-  Vichnevetsky,  Editor),  pp.  46-68. 


68 


The  Laboratory  for  Numerical  Analj^is  is  an  integral  part  of  the  Institute  for  I%ysical 
Science  and  Technology  of  the  University  of  Maryland,  under  the  general  administration  of  the 
Director,  Institute  for  Physical  Science  and  Technology.  It  has  the  following  goals: 

•  To  conduct  research  in  the  mathematical  theory  and  computational  implementation  of 
numerical  anafysis  and  related  topics,  with  enqihasis  on  the  numerical  treatment  of 
linear  and  nonlinear  differential  equations  and  proUems  in  linear  and  nonlinear  algeln*a. 

•  To  help  bridge  gaps  between  computational  directions  in  engineering,  physics,  etc.,  and 
those  in  the  mathematical  communi^. 

•  To  provide  a  limited  consulting  service  in  all  areas  of  numerical  mathematics  to  the 
University  as  a  whole,  and  also  to  government  agencies  and  industries  in  the  State  of 
Maryland  and  the  Washington  Metropolitan  area. 

•  To  assist  with  the  education  of  numerical  anafysts,  eqpedalfy  at  the  postdoctoral  level, 
in  conjunction  with  the  Interdisciplinary  Aj^lied  Mathematics  Program  and  the 
programs  of  the  Mathematics  and  Computer  Science  Departments.  This  includes  active 
collaboration  with  govenunent  agencies  such  as  the  National  Institute  of  Standards  and 
Technology. 

•  To  be  an  international  center  of  study  and  research  for  foreign  students  in  numerical 

mathematics  who  are  supported  1^  foreign  governments  or  exchange  agencies 
(Fulbright,  etc.).  j 
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Further  information  may  be  obtained  from  Pnrf’essor  L  BabuSka,  Chairman,  Laboratory  for 
Numerical  Anafysis,  Institute  for  Physical  Science  and  Technology,  University  of  Maryland,  College 
Park,  Maryland  20742-2431. 


